Re: [R] installation of packages

[Prof Brian Ripley]

Please see the discussion in the rw-FAQ.

On Wed, 15 Aug 2007, [EMAIL PROTECTED] wrote:

> Dear All,
>
> Have just installed v2.5.1 on Windows XP. Works fine but I had quite a few
> pakages loaded for 2.5.0 (from contributed) and was wondering how I can
> get 2.5.1 to recognise them without having to reinstall them all.
>
> Is this possible or do I have to reinstall all the packages again?
>
> I required 2.5.1 for lme4 and matrix.
>
> Many thanks in advance.
>
> 
>
> Regards
>
> Robin Dobos,
> Livestock Research Officer (Livestock Production Systems),
> Beef Industry Centre of Excellence,
> NSW Department of Primary Industries,
> Armidale, NSW, Australia, 2351
>
> ph:  +61 2 6770 1824
> fax:  +61 2 6770 1830
> mobile: 0431 391 885
> email: [EMAIL PROTECTED]
>
> If we knew what it was we were doing,
> it would not be called research, would it?
>
> Albert Einstein
>
>
> This message is intended for the addressee named and may con...{{dropped}}
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] installation of packages

[robin . dobos]

Dear All,

Have just installed v2.5.1 on Windows XP. Works fine but I had quite a few 
pakages loaded for 2.5.0 (from contributed) and was wondering how I can 
get 2.5.1 to recognise them without having to reinstall them all.

Is this possible or do I have to reinstall all the packages again?

I required 2.5.1 for lme4 and matrix.

Many thanks in advance.



Regards

Robin Dobos, 
Livestock Research Officer (Livestock Production Systems),
Beef Industry Centre of Excellence,
NSW Department of Primary Industries,
Armidale, NSW, Australia, 2351

ph:  +61 2 6770 1824
fax:  +61 2 6770 1830
mobile: 0431 391 885
email: [EMAIL PROTECTED]

If we knew what it was we were doing, 
it would not be called research, would it?

Albert Einstein


This message is intended for the addressee named and may con...{{dropped}}

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] OT - use of R

[Chung-hong Chan]

A quick search of google scholar found 2556 cites to the classic paper
by Ihaka et al.

http://scholar.google.com/scholar?hl=en&lr=&cites=15992947024900415641



On 8/15/07, Andy Bunn <[EMAIL PROTECTED]> wrote:
> Hello, is there an up-to-date reference for how many people use R? I'm
> giving an R demo and want to cite wonderful R usage stats. How many
> people use it (or download it)? How often is R used in peer-reviewed
> pubs, etc. Is there any whiz-bang citation that says something like "R
> is great and developed by the best minds in statistical computing." Any
> thoughts?
>
> -Andy
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>


-- 
"The scientists of today think deeply instead of clearly. One must be
sane to think clearly, but one can think deeply and be quite insane."
Nikola Tesla
http://www.macgrass.com

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Mann-Whitney U

[J Dougherty]

Natalie,

It's best to provide at least a sample of your data.  Your field names suggest 
that your data might be collected in units of mm^2 or some similar 
measurement of area.  Why do you want to use Mann-Whitney, which will rank 
your data and then use those ranks rather than your actual data?  Unless your 
sample is quite small, why not use a two sample t-test?  Also,are your 
samples paired?  If they aren't, did you use the "paired = FALSE" option?

JWDougherty

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] X11 problems

[De-Jian,ZHAO]

Hi Pau,

The error message indicates that the Font Path is perhaps wrong.
I think you should set the font path. I do not use Linux. I never 
encounter this problem in Windows XP.

About the warning messages, you can try "?locales" to get some
information. These warnings do not prevent the function from
working.

--
On Wed, Aug 15, 2007 10:03, Pau Marc Munoz Torres wrote:
> Hi
>
>  I'm working in a ubuntu feisty OS, when I try to start X11() i
get
> the
> following message
>
>> X11()
> Error in X11() : could not find any X11 fonts
> Check that the Font Path is correct.
> In addition: Warning messages:
> 1: locale not supported by Xlib: some X ops will operate in C
locale
> 2: X cannot set locale modifiers
>
> Can some body tell me what to do?
>
>
> --
> Pau Marc Mu�oz Torres
>
> Laboratori de Biologia Computacional
> Institut de  Biotecnologia   i Biomedicina Vicent
> Villar
> Universitat Autonoma de Barcelona
> E-08193 Bellaterra (Barcelona)
>
> tel�fon: 93 5812807
> Email : [EMAIL PROTECTED]
>
>   [[alternative HTML version deleted]]
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] X11 problems

[Marc Schwartz]

On Wed, 2007-08-15 at 04:03 +0200, Pau Marc Munoz Torres wrote:
> Hi
> 
>  I'm working in a ubuntu feisty OS, when I try to start X11() i get the
> following message
> 
> > X11()
> Error in X11() : could not find any X11 fonts
> Check that the Font Path is correct.
> In addition: Warning messages:
> 1: locale not supported by Xlib: some X ops will operate in C locale
> 2: X cannot set locale modifiers
> >
> 
> Can some body tell me what to do?

There seems to be a fair number of posts pertaining to issues with X11,
fonts and Ubuntu. As I use Fedora, there may be some subtleties in play
here that I am missing.

However, the message also suggests a problem with your locale, which has
also been referred to in posts and is referenced in the R Installation
and Administration Manual in Appendix C.1 "X11 Issues".

A quick workaround may be to start R, but set LANG=C on the command line
first:

  $ LANG=C R

Then try to start the X11() device and see if you still get the
messages.

If the locale message is gone, but the font path related messages are
still present, use:

  RSiteSearch("X11 Ubuntu")

which will present you with some possibilities for fixing the xorg.conf
file relative to defining font path settings. Seems that this is
somewhat common on Ubuntu upgrades where the font path changed from one
version to another, at least from Dapper to Edgy anyway.

HTH,

Marc Schwartz

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] vertically oriented color key in heatmaps

[adrian]

Rajarshi Guha <[EMAIL PROTECTED]> writes:

> Hi, I have some data which I was plotting using image(). I wanted to
> add a vertical color key to the plot and I found that heatmap.2 in
> gplots does let me add a color key. However, I was thinking of a
> vertical bar with the color range rather than  the style that gplots
> provides.

In package 'spatstat' this is the default behaviour for plotting a pixel
image (object of class "im").
If M is your matrix of data

 library(spatstat)
 Z <- as.im(t(M))
 plot(Z)
 plot(Z, col=topo.colors(128))

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] X11 problems

[Pau Marc Munoz Torres]

Hi

 I'm working in a ubuntu feisty OS, when I try to start X11() i get the
following message

> X11()
Error in X11() : could not find any X11 fonts
Check that the Font Path is correct.
In addition: Warning messages:
1: locale not supported by Xlib: some X ops will operate in C locale
2: X cannot set locale modifiers
>

Can some body tell me what to do?


-- 
Pau Marc Muñoz Torres

Laboratori de Biologia Computacional
Institut de  Biotecnologia   i Biomedicina Vicent
Villar
Universitat Autonoma de Barcelona
E-08193 Bellaterra (Barcelona)

telèfon: 93 5812807
Email : [EMAIL PROTECTED]

[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] simulate data from multivariate normal with pre-specified correlation matrix

[Thomas Harte]

# here's your example correlation matrix:
sigma<- matrix(c(1.00, 0.75, 0,
 0.75, 1.00, 0,
 0.00, 0.00, 0), nr=3, byrow=TRUE)
chol(sigma)
# Error in chol(sigma) : the leading minor of order 3 is not positive definite

# DUH!

# let's chop off that dangling row and column:
sigma<- sigma[1:2, 1:2]
N<- 1000
set.seed(1)
samp<- matrix(rnorm(N*nrow(sigma)), N, nrow(sigma))%*% chol(sigma)
# check if the empirical correlation is close to the theoretical sigma:
cor(samp)
#  [,1]  [,2]
#[1,] 1.000 0.7502607
#[2,] 0.7502607 1.000

# not bad!



>   Message: 50
>   Date: Mon, 13 Aug 2007 17:15:12 -0400
>   From: Jiao Yang <[EMAIL PROTECTED]>
>   Subject: [R] simulate data from multivariate normal with pre-specified
>   correlation matrix
>   To: r-help@stat.math.ethz.ch
>   Message-ID: <[EMAIL PROTECTED]>
>   Content-Type: text/plain; charset=us-ascii

>   For example, the correlation matrix is 3x3 and looks like 

>   1   0.75  0  0  0
>   0.751 0  0  0
>0   0 0   0  0

>   Can I write the code like  this?

>   p<- 3 >  number of variables per observation
>   N<- 10 >  number of samples

>   >  define population correlation matrix sigma
>   sigma<-matrix(0,p,p) > creates a px p matrix of 0
>   rank<-2
>   for (i in 1:rank){
>   for (j in 1:rank){
>   rho<-0.75
>sigma[i,j]<-rho^abs(i-j)
>   sigma[i,i]<-1
>}
>   }

>   >  create a data set of 10 observations from multivariate normalcovariance 
> matrix
sigma  
>   library(MASS)
>   Xprime<-mvrnorm(N,rep(1,p), sigma )

>   Also, if I calculate

>   S<- cov(Xprime) 

>   will S be the sample correlation matrix?

>   Thank you very much for your time!

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Mann-Whitney U

[Peter Dalgaard]

Prof Brian Ripley wrote:
> On Tue, 14 Aug 2007, Natalie O'Toole wrote:
>
>   
>> Hi,
>>
>> Could someone please tell me how to perform a Mann-Whitney U test on a
>> dataset with 2 groups where one group has more data values than another?
>>
>> I have split up my 2 groups into 2 columns in my .txt file i'm using with
>> R. Here is the code i have so far...
>>
>> group1 <- c(LeafArea2)
>> group2 <- c(LeafArea1)
>> wilcox.test(group1, group2)
>>
>> This code works for datasets with the same number of data values in each
>> column, but not when there is a different number of data values in one
>> column than another column of data.
>> 
>
> There is an example of that scenario on the help page for wilcox.test, so 
> it does 'work'.  What exactly went wrong for you?
>
>   
>> Is the solution that i have to have a null value in the data column with
>> the fewer data values?
>>
>> I'm testing for significant diferences between the 2 groups, and the
>> result i'm getting in R with the uneven values is different from what i'm
>> getting in SPSS.
>> 
>
> We need a worked example.  As the help page says, definitions do differ. 
> If you can provide a reproducible example in R and the output from SPSS we 
> may be able to tell you how to relate that to what you see in R.
>
> [...]
>
>   
>> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>> 
>
> As it says, we really need such code (and the output you get) to be able 
> to help you.
>
>   
Also, "two variables of different length in two columns" is not a good 
idea. If you read in things in parallel columns, it would usually imply 
paired data. If one column is shorter, you may be reading different data 
than you think. Check e.g. the "sleep" data for a better format.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Can I calculcate the percentage of a gamma function area below a cutoff value?

[Roland Rau]

Hi

Jake Verschuyl wrote:
> Hi there,
> 
>  
> 
> I have some bird flight height data that follows a gamma distribution.  The
> data (x) goes from 0 to 700 meters (n=1055).  The calculated parameters
> calculated from the fitdistr(x) are (shape = 5.1379, rate = 0.017541), and
> therefore the scale (1/rate) = 57.00929.  I would like to calculate the
> percentage of the function area that occurs below 50 meters (0-50m).  
> 

I guess
?pgamma
will help you?

and in your case:
pgamma(q=50, shape = 5.1379, rate = 0.017541)
 > pgamma(q=50, shape = 5.1379, rate = 0.017541)
[1] 0.001620588
 > pgamma(700, shape = 5.1379, rate = 0.017541)
[1] 0.9927367

Was this what you were looking for?

 >
 > P Please consider the environment before printing this email.
 >
 >
Don't worry, I will not print it. And in return:
Please consider reading the posting guide...
where it is stated:
=> Read at least the relevant section in "An Introduction to R"

There is a section (Section 8) which answers (almost) all questions on 
probability distributions and R.

Best,
Roland

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] OT - use of R

[Andy Bunn]

Hello, is there an up-to-date reference for how many people use R? I'm
giving an R demo and want to cite wonderful R usage stats. How many
people use it (or download it)? How often is R used in peer-reviewed
pubs, etc. Is there any whiz-bang citation that says something like "R
is great and developed by the best minds in statistical computing." Any
thoughts?

-Andy

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] CDF of a gamma distribution

[Rolf Turner]

On 15/08/2007, at 9:53 AM, Jake Verschuyl wrote:

> I have some bird flight height data that follows a gamma  
> distribution.  The
> data (x) goes from 0 to 700 meters (n=1055).  The calculated  
> parameters
> calculated from the fitdistr(x) are (shape = 5.1379, rate =  
> 0.017541), and
> therefore the scale (1/rate) = 57.00929.  I would like to calculate  
> the
> cumulative density (or proportion of the curve) of the function at  
> x= 50
> meters.
>
> Is that possible, and are there any functions that will help?

?pgamma

##
Attention:\ This e-mail message is privileged and confidenti...{{dropped}}

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] ANOVA: Factorial designs with a separate control

[lamack lamack]

Dear all, I would like to run in R the anova showed in the following 
pamphlet.

http://www.for.gov.bc.ca/hre/biopamph/pamp14.pdf

For A = 0 and B =0 I have de control group.

Best regards.

AB   replication   response
0 0 1 24
0 0 2 27
0 0 3 36
0 0 4 28
0 0 5 32
1 1 1 43
1 1 2 39
1 1 3 32
1 1 4 36
1 1 5 50
1 2 1 34
1 2 2 35
1 2 3 45
1 2 4 37
1 2 5 52
1 3 2 34
1 3 3 35
1 3 4 58
1 3 5 35
1 4 1 26
1 4 2 49
1 4 3 44
1 4 4 39
1 4 5 37
1 5 1 37
1 5 2 32
1 5 3 33
1 5 4 37
1 5 5 36
2 1 1 43
2 1 2 33
2 1 3 33
2 1 4 38
2 1 5 41
2 2 1 42
2 2 2 45
2 2 3 34
2 2 4 29
2 2 5 33
2 3 1 33
2 3 2 42
2 3 3 30
2 3 4 28
2 3 5 40
2 4 1 39
2 4 2 33
2 4 3 20
2 4 4 29
2 4 5 24
2 5 1 37
2 5 2 46
2 5 3 33
2 5 4 27
2 5 5 34

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] CDF of a gamma distribution

[Francisco J. Zagmutt]

pgamma(q=50,shape=5.1379,rate=0.017541)

Francisco




Jake Verschuyl wrote:
>  
> 
> Hi there,
> 
>  
> 
> I have some bird flight height data that follows a gamma distribution.  The
> data (x) goes from 0 to 700 meters (n=1055).  The calculated parameters
> calculated from the fitdistr(x) are (shape = 5.1379, rate = 0.017541), and
> therefore the scale (1/rate) = 57.00929.  I would like to calculate the
> cumulative density (or proportion of the curve) of the function at x= 50
> meters.  
> 
>  
> 
> Is that possible, and are there any functions that will help?  
> 
> Thanks!
> 
>  
> 
>  
> 
> Jake Verschuyl
> Ecologist
> 
> P.O. Box 2561, Mount Vernon, WA 98273
> 
> Tel: 360.899.5156
> 
> Fax: 360.899.5146 
> 
> P Please consider the environment before printing this email.
> 
>  
> 
> 
>   [[alternative HTML version deleted]]
> 
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Can I calculcate the percentage of a gamma function area below a cutoff value?

[Jake Verschuyl]

Hi there,

 

I have some bird flight height data that follows a gamma distribution.  The
data (x) goes from 0 to 700 meters (n=1055).  The calculated parameters
calculated from the fitdistr(x) are (shape = 5.1379, rate = 0.017541), and
therefore the scale (1/rate) = 57.00929.  I would like to calculate the
percentage of the function area that occurs below 50 meters (0-50m).  

 

Is that possible?  

Thanks!

 

 

Jake Verschuyl
Ecologist

P.O. Box 2561, Mount Vernon, WA 98273

Tel: 360.899.5156

Fax: 360.899.5146 

P Please consider the environment before printing this email.

 


[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] CDF of a gamma distribution

[Jake Verschuyl]

 

Hi there,

 

I have some bird flight height data that follows a gamma distribution.  The
data (x) goes from 0 to 700 meters (n=1055).  The calculated parameters
calculated from the fitdistr(x) are (shape = 5.1379, rate = 0.017541), and
therefore the scale (1/rate) = 57.00929.  I would like to calculate the
cumulative density (or proportion of the curve) of the function at x= 50
meters.  

 

Is that possible, and are there any functions that will help?  

Thanks!

 

 

Jake Verschuyl
Ecologist

P.O. Box 2561, Mount Vernon, WA 98273

Tel: 360.899.5156

Fax: 360.899.5146 

P Please consider the environment before printing this email.

 


[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Panel data and imputed datasets

[Nathan Paxton]

Hi all,

I am hardly an expert, so I expect that this code is not the easiest/ 
most efficient way of getting where I want. Any suggestions in that  
direction would also be helpful.

I am working on panel analysis with five imputed datasets, generated  
by Amelia. To do panel analysis, it seemed that the plm package was  
the best, providing a convenient wrapper for fixed and random effects  
and Hausman tests. And since I've got five datasets, it seemed  
reasonable that what I wanted to do was to combine these as a big  
imputation data frame.  The plm regression function requires a data  
argument, and simply using the imputationList result didn't work, nor  
did calling the names(allhiv) work (that resulted in "imputations"  
and "call").

Where do I point the plm command to look for the data? Or have I  
gone about this in entirely the wrong way?

Thanks in advance for the help.

Best,
-N

Here's the code I have so far:
#Read in the imputed datasets
data1 <- read.dta("outdata1.dta")
data2 <- read.dta("outdata2.dta")
data3 <- read.dta("outdata3.dta")
data4 <- read.dta("outdata4.dta")
data5 <- read.dta("outdata5.dta")

#Set the datasets as panels (for subsequent analysis in plm)
pdata.frame(data1, "country", time="year")
pdata.frame(data2, "country", time="year")
pdata.frame(data3, "country", time="year")
pdata.frame(data4, "country", time="year")
pdata.frame(data5, "country", time="year")

#Using mitools, combine the panel data frames 1-5
allhiv <- imputationList(list(data1,data2,data3,data4,data5))

#A test regression
form2<-APIPolSupport~subnatexpoftotal+adultprev.lag+pressfreedom

#The following is something I tried but which didn't work.
estimateHIV1 <-with(allhiv, plm(form2, data=imputations[1:5,]))

--
Nathan A. Paxton
Ph.D. Candidate
Dept. of Government, Harvard University

Resident Tutor
John Winthrop House, Harvard University

napaxton AT fas DOT harvard DOT edu
http://www.fas.harvard.edu/~napaxton
 
===
When you have to stay eight years away from California, you live in a  
perpetual state of homesickness.
 - Ronald Reagan

The most courageous act is still to think for yourself.  Aloud.
 -Coco Chanel

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Mann-Whitney U

[Prof Brian Ripley]

On Tue, 14 Aug 2007, Natalie O'Toole wrote:

> Hi,
>
> Could someone please tell me how to perform a Mann-Whitney U test on a
> dataset with 2 groups where one group has more data values than another?
>
> I have split up my 2 groups into 2 columns in my .txt file i'm using with
> R. Here is the code i have so far...
>
> group1 <- c(LeafArea2)
> group2 <- c(LeafArea1)
> wilcox.test(group1, group2)
>
> This code works for datasets with the same number of data values in each
> column, but not when there is a different number of data values in one
> column than another column of data.

There is an example of that scenario on the help page for wilcox.test, so 
it does 'work'.  What exactly went wrong for you?

> Is the solution that i have to have a null value in the data column with
> the fewer data values?
>
> I'm testing for significant diferences between the 2 groups, and the
> result i'm getting in R with the uneven values is different from what i'm
> getting in SPSS.

We need a worked example.  As the help page says, definitions do differ. 
If you can provide a reproducible example in R and the output from SPSS we 
may be able to tell you how to relate that to what you see in R.

[...]

> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

As it says, we really need such code (and the output you get) to be able 
to help you.

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] GML with tweedie: AIC=NA

[Gordon Smyth]

Dear Catarina,

As I've already indicated, I have no plans to include AIC for tweedie 
families, because the AIC computation is much more computationally 
intensive and less reliable than the actually fitting of the glm. 
That fact that dtweedie() is not currently available for Windows is a 
demonstration of my point.

I'm sorry that this causes you a problem. I don't pretend to know 
your problem but, in my experience, it can often be an advantage to 
do your model selection more hands on rather than relying on an 
automatic algorithm.

Best wishes
Gordon

At 09:06 PM 14/08/2007, Catarina Miranda wrote:
>Thanks for your reply, but I still couldn't solve the problem.
>
>I am using the package statmod, and I need the AIC because I want to 
>use the step function (I am modelling many species, so I would 
>prefer to do the step automatically).
>I can't find the tweedie  package in the R packages list and I don't 
>know how to download it from the package source available in 
>http://cran.r-project.org/src/contrib/Descriptions/tweedie.html
> 
>.
>However, after running the  tweedie.R file from that package source 
>I am able to use the dtweedie() function, but still I didn't figure 
>out a way to do the step (or to get the AIC from the glm command).
>
>Thank you again for your help;
>
>Catarina
>
>
>On 13/08/07, Gordon Smyth 
><[EMAIL PROTECTED] > wrote:
>Dear Catarina,
>
>I prefer to leave the AIC value as NA for the tweedie GLM family
>because it takes extra time to compute and is only occasionally
>wanted. It's easy to compute the AIC yourself using the dtweedie()
>function of the tweedie package.
>
>Best wishes
>Gordon
>
>At 03:05 AM 14/08/2007, Catarina Miranda wrote:
> >Dear Gordon;
> >
> >I have also sent this email to R help mailing list, so I apologize
> >for duplicated mailing.
> >I am modelling densities of some species of birds, and I have a
> >problem with a great amount of zeros.
> >I have decided to try GLMs with the tweedie family, but in all the
> >models I have tried  I got an NA for the AIC value.
> >Just  to check the problem I've compared the a glm using the
> >Gaussian family with the identity link and a glm using the tweedie
> >family with var.power=0 and link.power=1. These are equal, as
> >expected, except the fact that the tweedie output gives me an NA 
> for the AIC.
> >Could you help me with this problem?
> >Below you can find the two outputs I refer.
> >
> >Best Wishes;
> >
> >Catarina
> >
> > > summary(glm(formula=ACIN~DIST_REF+DIST_H2O+DIST_OST+
> > COTA+H2O_SUP+vasa,family=gaussian(link="identity")))
> >Call:glm(formula = ACIN ~ DIST_REF + DIST_H2O + DIST_OST + COTA
> >+ H2O_SUP + vasa, family = gaussian(link = "identity"))
> >Deviance
> >Residuals:   Min 1Q Median 3QMax
> >-0.112792   -0.042860  -0.021113  -0.006311   1.551824
> >Coefficients:  Estimate Std. Error t value
> >Pr(>|t|)  (Intercept)
> >-6.625e-02  5.454e-02  -1.215   0.2256  DIST_REF 3.581e-06
> >1.336e-05   0.268   0.7889  DIST_H2O-
> > 3.168e-05  1.527e-05  -2.074   0.0391
> >*DIST_OST-1.799e-05  1.953e-05  -0.921   0.3579  COTA
> >5.648e-04  2.470e-04   2.287   0.0230
> >*H2O_SUP -2.172e-04  3.994e-04  -0.544   0.5870  vasa
> >3.695e-02   4.573e-020.808   0.4199  ---Signif. codes:  0 '***'
> >0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
> >(Dispersion parameter for gaussian family taken to be 0.02151985)
> > Null deviance: 5.6028   on 257  degrees of freedomResidual
> > deviance: 5.4015  on 251  degrees of freedomAIC: -249.33
> >Number of Fisher Scoring iterations: 2
> >
> >
> > > summary(glm(formula=ACIN~DIST_REF+DIST_H2O+DIST_OST+
> > COTA+H2O_SUP+vasa,control=
> > glm.control(maxit=750),family=tweedie(var.power=0, link.power=1)))
> >Call:glm(formula = ACIN ~ DIST_REF + DIST_H2O + DIST_OST + COTA
> >+ H2O_SUP + vasa, family = tweedie( var.power = 0, link.power =
> >1), control = glm.control (maxit = 750))
> >Deviance
> >Residuals:   Min 1Q Median 3QMax
> >-0.112792  -0.042860  -0.021113  -0.006311   1.551824
> >Coefficients:  Estimate Std. Error t value
> >Pr(>|t|)  (Intercept) -
> >6.625e-02  5.454e-02  -1.215   0.2256  DIST_REF 3.581e-06
> >1.336e-05   0.268   0.7889  DIST_H2O-3.168e-05   1.527e-05
> >-2.074   0.0391
> >*DIST_OST-1.799e-05  1.953e-05  -0.921   0.3579  COTA
> >5.648e-04  2.470e-042.287   0.0230
> >*H2O_SUP -2.172e-04  3.994e-04  -0.544   0.5870  vasa
> >3.695e-02   4.573e-02   0.808   0.4199  ---Signif. codes:  0 '***'
> >0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
> >(Dispersion parameter for Tweedie family taken to be 0.02151985)
> > Null deviance: 5.6028  on 257  degrees of freedomResidual
> > deviance: 5.4015  on 251  degrees of freedomAIC: NA
> >Number of Fisher Scoring iterations: 2
> >
> >
>

__
R-help@stat.math.ethz.

Re: [R] Mann-Whitney U

[Marc Schwartz]

On Tue, 2007-08-14 at 14:45 -0600, Natalie O'Toole wrote:
> Hi,
> 
> Could someone please tell me how to perform a Mann-Whitney U test on a 
> dataset with 2 groups where one group has more data values than another?
> 
> I have split up my 2 groups into 2 columns in my .txt file i'm using with 
> R. Here is the code i have so far...
> 
> group1 <- c(LeafArea2)
> group2 <- c(LeafArea1)
> wilcox.test(group1, group2)
> 
> This code works for datasets with the same number of data values in each 
> column, but not when there is a different number of data values in one 
> column than another column of data.
> 
> Is the solution that i have to have a null value in the data column with 
> the fewer data values?
> 
> I'm testing for significant diferences between the 2 groups, and the 
> result i'm getting in R with the uneven values is different from what i'm 
> getting in SPSS.
> 
> Help please!
> 
> Nat

You will need to provide any error messages that you are getting. There
is a two sample example in ?wilcox.test that shows that the function can
handle two vectors with differing sizes.

Having the output of str(group1) and str(group2) may also prove useful.

You may also wish to pay attention to the "Note" in ?wilcox.test which,
if you are getting differing results between SPSS and R, may provide
some insight into why, presuming that you can gain the same information
about SPSS.

HTH,

Marc Schwartz

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Convert factor to numeric vector of labels

[Marc Schwartz]

I think that you grossly underestimate the frequency of use of factors
in R, not to mention that factors are stored more efficiently than
character vectors.

All modeling functions depend upon them.  Most testing, grouping and
plotting functions (base R and Lattice) either use them directly as
arguments or coerce character vectors to factors internally.

So, no...I would not advocate modifying such fundamental behavior.

UseRs should read the documentation before "jumping in with both feet"
so that they understand the underlying design philosophy behind R and
the actual documented functional behaviors. This would be superior to
moving forward with functional expectations that are predicated on false
assumptions and importantly, save you time.

In Falk's case, it seems reasonable, without having seen any actual
data, that the presumptive numeric column that was converted to a
factor, had non-numeric characters in it. 

Thus, that a numeric column was coerced to a factor on import should
have raised a red flag pointing to a data quality problem.  

Had the default behavior been otherwise, it is likely that Falk would
have proceeded with subsequent analyses without being aware of this
issue, perhaps resulting in a bad outcome.

The function that handles this in the read.table() family of functions
is called type.convert(). An example may be helpful:

Vec <- as.character(1:10)

> Vec
 [1] "1"  "2"  "3"  "4"  "5"  "6"  "7"  "8"  "9"  "10"

# Default behavior converts the character vector
# to numeric
> str(type.convert(Vec))
 int [1:10] 1 2 3 4 5 6 7 8 9 10


# Now add in a non-numeric character (ie. bad data)
Vec1 <- c(Vec, "a")

> Vec1
 [1] "1"  "2"  "3"  "4"  "5"  "6"  "7"  "8"  "9"  "10" "a" 


> str(type.convert(Vec1))
 Factor w/ 11 levels "1","10","2","3",..: 1 3 4 5 6 7 8 9 10 2 ...

Voilà

HTH,

Marc Schwartz


On Tue, 2007-08-14 at 13:47 -0600, Matthew Keller wrote:
> Hi all,
> 
> If we, the R community, are endeavoring to make R user friendly
> (gasp!), I think that one of the first places to start would be in
> setting stringsAsFactors = FALSE. Several times I've run into
> instances of folks decrying R's "rediculous usage of memory" in
> reading data, only to come to find out that these folks were
> unknowingly importing certain columns as factors. The fix is easy once
> you know it, but it isn't obvious to new users, and I'd bet that it
> turns some % of people off of the program. Factors are not used often
> enough to justify this default behavior in my opinion. When factors
> are used, the user knows to treat the variable as a factor, and so it
> can be done on a case-by-case (or should I say variable-by-variable?)
> basis.
> 
> Is this a default that should be changed?
> 
> Matt
> 
> 
> On 8/13/07, John Kane <[EMAIL PROTECTED]> wrote:
> > This is one of R's rather _endearing_  little
> > idiosyncrasies. I ran into it a while ago.
> > http://finzi.psych.upenn.edu/R/Rhelp02a/archive/98090.html
> >
> >
> > For some reason, possibly historical, the option
> > "stringAsFactors" is set to TRUE.
> >
> > As Prof Ripley says FAQ 7.10 will tell you
> > as.numeric(as.character(f)) # for a one-off conversion
> >
> > >From Gabor Grothendieck  A one-off solution for a
> > complete data.frame
> >
> > DF <- data.frame(let = letters[1:3], num = 1:3,
> >  stringsAsFactors = FALSE)
> >
> > str(DF)  # to see what has happened.
> >
> > You can reset the option globally, see below.  However
> > you might want to read Gabor Grothendieck's comment
> > about this in the thread referenced above since it
> > could cause problems if you transfer files alot.
> >
> > Personally I went with the global option since I don't
> > tend to transfer programs to other people and I was
> > getting tired of tracking down errors in my programs
> > caused by numeric and character variables suddenly
> > deciding to become factors.
> >
> > >From Steven Tucker:
> >
> > You can also this option globally with
> >  options(stringsAsFactors = TRUE)  # in
> > \library\base\R\Rprofile
> >
> > --- Falk Lieder <[EMAIL PROTECTED]> wrote:
> >
> > > Hi,
> > >
> > > I have imported a data file to R. Unfortunately R
> > > has interpreted some
> > > numeric variables as factors. Therefore I want to
> > > reconvert these to numeric
> > > vectors whose values are the factor levels' labels.
> > > I tried
> > > as.numeric(),
> > > but it returns a vector of factor levels (i.e.
> > > 1,2,3,...) instead of labels
> > > (i.e. 0.71, 1.34, 2.61,…).
> > > What can I do instead?
> > >
> > > Best wishes, Falk
> >

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Mann-Whitney U

[Natalie O'Toole]

Hi,

Could someone please tell me how to perform a Mann-Whitney U test on a 
dataset with 2 groups where one group has more data values than another?

I have split up my 2 groups into 2 columns in my .txt file i'm using with 
R. Here is the code i have so far...

group1 <- c(LeafArea2)
group2 <- c(LeafArea1)
wilcox.test(group1, group2)

This code works for datasets with the same number of data values in each 
column, but not when there is a different number of data values in one 
column than another column of data.

Is the solution that i have to have a null value in the data column with 
the fewer data values?

I'm testing for significant diferences between the 2 groups, and the 
result i'm getting in R with the uneven values is different from what i'm 
getting in SPSS.

Help please!

Nat



 

This communication is intended for the use of the recipient to which it is 
addressed, and may
contain confidential, personal, and or privileged information. Please 
contact the sender
immediately if you are not the intended recipient of this communication, 
and do not copy,
distribute, or take action relying on it. Any communication received in 
error, or subsequent
reply, should be deleted or destroyed.

 

This communication is intended for the use of the recipient to which it is 
addressed, and may
contain confidential, personal, and or privileged information. Please 
contact the sender
immediately if you are not the intended recipient of this communication, 
and do not copy,
distribute, or take action relying on it. Any communication received in 
error, or subsequent
reply, should be deleted or destroyed.
[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] State Space Modelling

[sj]

If your looking at state space models for forecasting, the Hyndman's
forecasting has quite a bit.

On 8/14/07, Bernardo Ribeiro <[EMAIL PROTECTED]> wrote:
>
> Here is one of the codes.
>
> Thanks a lot...
>
>
> module(finmetrics)
>
> genvasicek.ssf = function(param, tau=NULL, freq=1/12){
>
> ## 1. Check for valid inputs
>
> if(length(param) < 5)
> stop("Parameters must have length greater than 4")
> N = length(param) - 13
> if (length(tau) != N)
> stop("Length of Tau is inconsistent with Parameters")
>
> ## 2.Extract Parameters and Impose Constraints
>
> k1 = (param[1])
> k2 = (param[2])
> k3 = (param[3])
>
> lambda1 = param[4]
> lambda2 = param[5]
> lambda3 = param[6]
>
> Delta = (param[7])
>
> s1 = (param[8])
> s2 = (param[9])
> s3 = (param[10])
>
> r21 = param[11]
> r31 = param[12]
> r32 = param[13]
>
> st1 = (param[14])
> st2 = (param[15])
> st3 = (param[16])
> st4 = (param[17])
>
> ## 3. Compute Matrixes
>
> A1 = function(x,tau){
> (1-exp(-x*tau))/x
> }
>
> u1 = -A1(k1,tau)  ## ok
> u2 = -A1(k2,tau)  ## ok
> u3 = -A1(k3,tau)  ## ok
> H = cbind (-u1/tau,-u2/tau,-u3/tau)  ## ok
>
> v1 = (lambda1/k1)*(tau - A1(k1,tau)) + (lambda2/k2)*(tau - A1(k2,tau))
> +
> (lambda3/k3)*(tau - A1(k3,tau))
> v2 = +3 * Delta * tau  # Teste
>
> v311 = ((s1*s1)/(k1*k1)) * (tau - A1(k1,tau) - A1(k1,tau) +
> A1(2*k1,tau))
> v312 = ((s1*s2*r21)/(k1*k2)) * (tau - A1(k1,tau) - A1(k2,tau) +
> A1(k1+k2,tau))
> v313 = ((s1*s3*r31)/(k1*k3)) * (tau - A1(k1,tau) - A1(k3,tau) +
> A1(k1+k3,tau))
>
> v322 = ((s2*s2)/(k2*k2)) * (tau - A1(k2,tau) - A1(k2,tau) +
> A1(2*k2,tau))
> v321 = ((s1*s2*r21)/(k1*k2)) * (tau - A1(k1,tau) - A1(k2,tau) +
> A1(k1+k2,tau))
> v323 = ((s2*s3*r32)/(k2*k3)) * (tau - A1(k2,tau) - A1(k3,tau) +
> A1(k2+k3,tau))
>
> v333 = ((s3*s3)/(k3*k3)) * (tau - A1(k3,tau) - A1(k3,tau) +
> A1(2*k3,tau))
> v331 = ((s1*s3*r31)/(k1*k3)) * (tau - A1(k3,tau) - A1(k1,tau) +
> A1(k3+k1,tau))
> v323 = ((s2*s3*r32)/(k2*k3)) * (tau - A1(k2,tau) - A1(k3,tau) +
> A1(k2+k3,tau))
>
> v = v1+v2+0.5*(v311+v312+v313+v322+v321+v323+v333+v331+v323)
>
> d = -v
>
> A = matrix(c(0,0,0,0,0,0,0,0,0),ncol=3)
> A[1,1] = 1-k1*freq
> A[1,2] = 0
> A[1,3] = 0
> A[2,1] = 0
> A[2,2] = 1-k2*freq
> A[2,3] = 0
> A[3,1] = 0
> A[3,2] = 0
> A[3,3] = 1-k3*freq
>
> ct = c(-lambda1*freq,-lambda2*freq,-lambda3*freq)
>
> Q = matrix(c(0,0,0,0,0,0,0,0,0),ncol=3)
> Q[1,1] = s1^2
> Q[1,2] = s1*s2*r21
> Q[1,3] = s1*s3*r31
> Q[2,1] = s2*s1*r21
> Q[2,2] = s2^2
> Q[2,3] = s2*s3*r32
> Q[3,1] = s1*s3*r31
> Q[3,2] = s2*s3*r32
> Q[3,3] = s3^2
>
> ## 4. Compute the State Space Form
>
> mDelta = matrix(c(ct,d),ncol=1)
>
> mPhi = rbind(A,H)
>
> mOmega1 = matrix(c(Q,0,0,0,0,0,0,0,0,0,0,0,0),ncol=7)
> mOmega2 =
> matrix(c(0,0,0,0,0,0,0,0,0,0,0,0,diag(c(st1,st2,st3,st4))),ncol=7)
>
> mOmega = rbind(mOmega1,mOmega2)
>
> mSigma = matrix(c(0,0,0,0,0,0,0,0,0,0,0,0),ncol=3)
>
> mSigma[1,1] = ((s1^2)/(2*k1))*(1-exp(-2*k1*freq))
> mSigma[2,2] = ((s2^2)/(2*k2))*(1-exp(-2*k2*freq))
> mSigma[3,3] = ((s3^2)/(2*k3))*(1-exp(-2*k3*freq))
>
> ## 7. Return State Space Form
>
> ssf.mod = list(mDelta=mDelta, mPhi=mPhi, mOmega=mOmega)
> CheckSsf(ssf.mod)
>
> }
>
> start = c(0.1,0.1,0.1,0.3,0.3,0.3,0.06,0.005,0.005,0.005,0.001,0.001,0.001
> ,
> 0.003,0.1,0.003,0.01)
>
> names(start)=c("logkappa1","logkappa2","logkappa3","lambda1","lambda2","lambda3","logDelta","logSigma1","logSigma2","logSigma3","Corr12","Corr13","Corr23","logSigmat1","logSigmat2","logSigmat3","logSigmat4")
>
> tau = c(0.25,0.5,1,5)
>
> l = c(0,0,0,0,0,0,0,0,0,-0,-1,-1,-1,0,0,0,0)
> u = c(10,10,10,10,10,10,0.5,10,10,10,1,1,1,10,10,10,10)
>
>
> ans.vasicek = SsfFit(start, fama.bliss, genvasicek.ssf, tau=tau,
> freq=1/12,
> trace=T, lower=l, upper=u, control=nlminb.control(abs.tol=1e-6,,
> rel.tol=1e-6, x.tol=1e-6, eval.max=1000,iter.max=135))
>
>
> On 8/13/07, Bernardo Ribeiro <[EMAIL PROTECTED]> wrote:
> >
> > Hey all,
> >
> > I am trying to work under a State Space form, but I didn't get the help
> > exactly.
> > Have anyone eles used this functions?
> >
> > I was used to work with S-PLUS, but I have some codes I need to adpt.
> >
> > Thanks alot,
> >
> > Bernardo
> >
>
> [[alternative HTML version deleted]]
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R

Re: [R] Convert factor to numeric vector of labels

[Bert Gunter]

Matt:

I believe you have confused issues.

Setting stringsAsFactors = FALSE would dramatically **increase** the amount
of memory used for storing character vectors, which is what factors are for.
So your proposed solution does exactly the opposite of what you want.

The issue you are worried about is when numeric fields are somehow
interpreted as non-numeric. This can happen for a variety of reasons (stray
characters in numeric fields,quotes around numbers,...). The solution is not
to set a global default that does the opposite of what you want in its
intended use, but to read the documentation and either set the appropriate
arguments (perhaps colClasses of read.table) or fix the original data before
R reads it (e.g. remove quotes and stray characters). Failing that, the
"one-off" solutions given are the correct way to handle what is a data
problem, not an R problem.

However, I should add that there are arguments for making stringsAsFactors =
FALSE; search the archives for discussions why. The memory penalty will have
to be paid, of course.


Bert Gunter
Genentech Nonclinical Statistics


-Original Message-
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] On Behalf Of Matthew Keller
Sent: Tuesday, August 14, 2007 12:48 PM
To: John Kane
Cc: Falk Lieder; r-help@stat.math.ethz.ch
Subject: Re: [R] Convert factor to numeric vector of labels

Hi all,

If we, the R community, are endeavoring to make R user friendly
(gasp!), I think that one of the first places to start would be in
setting stringsAsFactors = FALSE. Several times I've run into
instances of folks decrying R's "rediculous usage of memory" in
reading data, only to come to find out that these folks were
unknowingly importing certain columns as factors. The fix is easy once
you know it, but it isn't obvious to new users, and I'd bet that it
turns some % of people off of the program. Factors are not used often
enough to justify this default behavior in my opinion. When factors
are used, the user knows to treat the variable as a factor, and so it
can be done on a case-by-case (or should I say variable-by-variable?)
basis.

Is this a default that should be changed?

Matt


On 8/13/07, John Kane <[EMAIL PROTECTED]> wrote:
> This is one of R's rather _endearing_  little
> idiosyncrasies. I ran into it a while ago.
> http://finzi.psych.upenn.edu/R/Rhelp02a/archive/98090.html
>
>
> For some reason, possibly historical, the option
> "stringAsFactors" is set to TRUE.
>
> As Prof Ripley says FAQ 7.10 will tell you
> as.numeric(as.character(f)) # for a one-off conversion
>
> >From Gabor Grothendieck  A one-off solution for a
> complete data.frame
>
> DF <- data.frame(let = letters[1:3], num = 1:3,
>  stringsAsFactors = FALSE)
>
> str(DF)  # to see what has happened.
>
> You can reset the option globally, see below.  However
> you might want to read Gabor Grothendieck's comment
> about this in the thread referenced above since it
> could cause problems if you transfer files alot.
>
> Personally I went with the global option since I don't
> tend to transfer programs to other people and I was
> getting tired of tracking down errors in my programs
> caused by numeric and character variables suddenly
> deciding to become factors.
>
> >From Steven Tucker:
>
> You can also this option globally with
>  options(stringsAsFactors = TRUE)  # in
> \library\base\R\Rprofile
>
> --- Falk Lieder <[EMAIL PROTECTED]> wrote:
>
> > Hi,
> >
> > I have imported a data file to R. Unfortunately R
> > has interpreted some
> > numeric variables as factors. Therefore I want to
> > reconvert these to numeric
> > vectors whose values are the factor levels' labels.
> > I tried
> > as.numeric(),
> > but it returns a vector of factor levels (i.e.
> > 1,2,3,...) instead of labels
> > (i.e. 0.71, 1.34, 2.61,.).
> > What can I do instead?
> >
> > Best wishes, Falk
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>


-- 
Matthew C Keller
Postdoctoral Fellow
Virginia Institute for Psychiatric and Behavioral Genetics

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] File Exchange?

[Dirk Eddelbuettel]

On Tue, Aug 14, 2007 at 02:34:45PM -0400, Wojciech Gryc wrote:
> Hi,
> 
> I've been using R for a number of years, and must say I really like it. This
> past summer, however, I started using MatLab for a course at school and
> found the file exchange it has (
> http://www.mathworks.com/matlabcentral/fileexchange/loadCategory.do) to be
> extremely useful.
> 
> I know that R has packages that people share, and that there's also RForge
> for more advanced tools... But is there anywhere where people share R
> scripts? For example, if I make a function that does one specific task, but
> I don't want to make a formal package, is there a centralized location where
> I can upload it and share it with others?

Try the R Wiki.

Dirk

-- 
Three out of two people have difficulties with fractions.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Convert factor to numeric vector of labels

[Matthew Keller]

Hi all,

If we, the R community, are endeavoring to make R user friendly
(gasp!), I think that one of the first places to start would be in
setting stringsAsFactors = FALSE. Several times I've run into
instances of folks decrying R's "rediculous usage of memory" in
reading data, only to come to find out that these folks were
unknowingly importing certain columns as factors. The fix is easy once
you know it, but it isn't obvious to new users, and I'd bet that it
turns some % of people off of the program. Factors are not used often
enough to justify this default behavior in my opinion. When factors
are used, the user knows to treat the variable as a factor, and so it
can be done on a case-by-case (or should I say variable-by-variable?)
basis.

Is this a default that should be changed?

Matt


On 8/13/07, John Kane <[EMAIL PROTECTED]> wrote:
> This is one of R's rather _endearing_  little
> idiosyncrasies. I ran into it a while ago.
> http://finzi.psych.upenn.edu/R/Rhelp02a/archive/98090.html
>
>
> For some reason, possibly historical, the option
> "stringAsFactors" is set to TRUE.
>
> As Prof Ripley says FAQ 7.10 will tell you
> as.numeric(as.character(f)) # for a one-off conversion
>
> >From Gabor Grothendieck  A one-off solution for a
> complete data.frame
>
> DF <- data.frame(let = letters[1:3], num = 1:3,
>  stringsAsFactors = FALSE)
>
> str(DF)  # to see what has happened.
>
> You can reset the option globally, see below.  However
> you might want to read Gabor Grothendieck's comment
> about this in the thread referenced above since it
> could cause problems if you transfer files alot.
>
> Personally I went with the global option since I don't
> tend to transfer programs to other people and I was
> getting tired of tracking down errors in my programs
> caused by numeric and character variables suddenly
> deciding to become factors.
>
> >From Steven Tucker:
>
> You can also this option globally with
>  options(stringsAsFactors = TRUE)  # in
> \library\base\R\Rprofile
>
> --- Falk Lieder <[EMAIL PROTECTED]> wrote:
>
> > Hi,
> >
> > I have imported a data file to R. Unfortunately R
> > has interpreted some
> > numeric variables as factors. Therefore I want to
> > reconvert these to numeric
> > vectors whose values are the factor levels' labels.
> > I tried
> > as.numeric(),
> > but it returns a vector of factor levels (i.e.
> > 1,2,3,...) instead of labels
> > (i.e. 0.71, 1.34, 2.61,…).
> > What can I do instead?
> >
> > Best wishes, Falk
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>


-- 
Matthew C Keller
Postdoctoral Fellow
Virginia Institute for Psychiatric and Behavioral Genetics

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] factorial Design with a separate control group

[lamack lamack]

Dear all, in R, how can I run a Factorial designs with a separate control?
Best regards.


#A = 0 and B = 0 == then control

ABrepresponse
0 0 1 24
0 0 2 27
0 0 3 36
0 0 4 28
0 0 5 32
1 1 1 43
1 1 2 39
1 1 3 32
1 1 4 36
1 1 5 50
1 2 1 34
1 2 2 35
1 2 3 45
1 2 4 37
1 2 5 52
1 3 2 34
1 3 3 35
1 3 4 58
1 3 5 35
1 4 1 26
1 4 2 49
1 4 3 44
1 4 4 39
1 4 5 37
1 5 1 37
1 5 2 32
1 5 3 33
1 5 4 37
1 5 5 36
2 1 1 43
2 1 2 33
2 1 3 33
2 1 4 38
2 1 5 41
2 2 1 42
2 2 2 45
2 2 3 34
2 2 4 29
2 2 5 33
2 3 1 33
2 3 2 42
2 3 3 30
2 3 4 28
2 3 5 40
2 4 1 39
2 4 2 33
2 4 3 20
2 4 4 29
2 4 5 24
2 5 1 37
2 5 2 46
2 5 3 33
2 5 4 27
2 5 5 34

_
Mande torpedos SMS do seu messenger para o celular dos seus amigos

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] [R-pkgs] 'fda' 1.2.2 is now available on CRAN.

[hadley wickham]

The fda package supports "Functional Data Analysis" and "Applied
Functional Data Analysis" by Bernard Silverman and James Ramsay.
Functional data analysis, which lots of us like to call "FDA", is
about the analysis of information on curves or functions. FDA is a
collection statistical techniques for answering questions like, "What
are the main ways in which the curves vary from one to another?" In
fact, most of the questions and problems associated with multivariate
data (PCA, LDA, clustering, ...) have functional counterparts. More
information about FDA can be found at
http://www.psych.mcgill.ca/misc/fda/.

This version (and the previous 1.2.1) includes bug fixes plus a
"scripts" subdirectory with R code to reproduce some of the analyses
in the two functional data analysis books by Ramsay and Silverman and
a "Continuously Stirred Tank Reactor (CSTR)" simulation discussed in a
Ramsay, et al., discussion paper to appear soon in the Journal of the
Royal Statistical Society-series B.

It also includes the draft of a presentation on "fda in Matlab & R"
(in PowerPoint and Adobe Acrobat PDF formats) for the UseR! 2007
conference this Friday, Aug. 10, 1:55 - 2:20 PM in Ames, IA.

Regards

Hadley Wickham
James Ramsey
Spencer Graves

___
R-packages mailing list
[EMAIL PROTECTED]
https://stat.ethz.ch/mailman/listinfo/r-packages

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] memory allocation glitches

[Peter Dalgaard]

Ben Bolker wrote:
> (not sure whether this is better for R-devel or R-help ...)
>
>   
Hardcore debugging is usually better off in R-devel. I'm leaving it in 
R-help though.

>  I am currently trying to debug someone else's package (they're
> not available at the moment, and I would like it to work *now*),
> which among other things allocates memory for a persistent
> buffer that gets used by various functions.
>
>  The first symptoms of a problem were that some things just
> didn't work under Windows but were (apparently) fine on Linux.
> I don't have all the development tools installed for Windows, so
> I started messing around under Linux, adding Rprintf() statements
> to the main code.
>
>  Once I did that, strange pointer-error-like inconsistencies started
> appearing -- e.g., the properties of some of the persistent variables
> would change if I did debug(function).  I'm wondering if anyone
> has any tips on how to tackle this -- figure out how to use valgrind?
> Do straight source-level debugging (R -d gdb etc.) and look for
> obvious problems?  The package uses malloc/realloc rather than
> Calloc/Realloc -- does it make sense to go through the code
> replacing these all and see if that fixes the problem?
>   
Valgrind is a good idea to try and as I recall it, the basic 
incantations are not too hard to work out (now exactly where is it that 
we wrote them down?). It only catches certain error types though, mostly 
use of uninitialized data and read/write off the ends of allocated 
blocks of memory.

If that doesn't catch it, you get to play with R -d  gdb. However, my 
experience is that line-by-line tracing is usually a dead end, unless 
you have the trouble spot pretty well narrowed down.

Apart from that, my usual procedure would be

1) find a minimal script reproducing the issue and hang onto it. Or at 
least as small as you can get it without losing the bug. Notice that any 
change to either the script or R itself may allow the bug to run away 
and hide somewhere else.

2) if memory corruption is involved, run under gdb, set a hardware 
watchpoint on the relevant location (this gets a little tricky sometimes 
because it might be outside the initial address space, in which case you 
need to somehow run the code for a while, break to gdb, and then set the 
watchpoint).

3) It is not unlikely that the watchpoint triggers several thousand 
times before the relevant one. You can conditionalize it; a nice trick 
is to use the gc_count.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] weights in lmer

[Sundar Dorai-Raj]

Try

weights = as.numeric(total)

BTW, there is a SIG (Special Interest Group) for lmer.

https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

HTH,

--sundar

Chris O'Brien said the following on 8/14/2007 11:00 AM:
> Dear R users,
> 
> Prof. Ripley just corrected my understanding of the use of weights in glm, 
> which I thought would allow me to correctly use lmer.  However I'm still 
> having problems.
> 
> My data takes the form of # of infected and uninfected individuals that 
> were measured over time under different treatments.  I'm using lmer to 
> adjust for the repeated measures over time.
> 
> In fitting the model:
> 
>  > model1=lmer(y~treatment+(time | ID),family=binomial,weights = total)
> 
> where y = proportion of animals infected (number infected/total)
>  total = number of infected + number uninfected
> 
> this returns the error:
> 
> Error in lmer(y ~ treatment + (time | ID), family = binomial, weights = 
> total) :
>  object `weights' of incorrect type
> 
> Can anyone tell me what's wrong with the weights?
> 
> cheers,
> Chris
> 
> 
> 
> 
> 
> Chris O'Brien
> Sonoran Desert Research Station and
> School of Natural Resources
> Biological Sciences, room 125
> University of Arizona
> Tucson, AZ 85721
> work (520) 623-3720
> [EMAIL PROTECTED]
> 
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] File Exchange?

[Wojciech Gryc]

Hi,

I've been using R for a number of years, and must say I really like it. This
past summer, however, I started using MatLab for a course at school and
found the file exchange it has (
http://www.mathworks.com/matlabcentral/fileexchange/loadCategory.do) to be
extremely useful.

I know that R has packages that people share, and that there's also RForge
for more advanced tools... But is there anywhere where people share R
scripts? For example, if I make a function that does one specific task, but
I don't want to make a formal package, is there a centralized location where
I can upload it and share it with others?

Thank you,
Wojciech

[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] memory allocation glitches

[Ben Bolker]

(not sure whether this is better for R-devel or R-help ...)

 I am currently trying to debug someone else's package (they're
not available at the moment, and I would like it to work *now*),
which among other things allocates memory for a persistent
buffer that gets used by various functions.

 The first symptoms of a problem were that some things just
didn't work under Windows but were (apparently) fine on Linux.
I don't have all the development tools installed for Windows, so
I started messing around under Linux, adding Rprintf() statements
to the main code.

 Once I did that, strange pointer-error-like inconsistencies started
appearing -- e.g., the properties of some of the persistent variables
would change if I did debug(function).  I'm wondering if anyone
has any tips on how to tackle this -- figure out how to use valgrind?
Do straight source-level debugging (R -d gdb etc.) and look for
obvious problems?  The package uses malloc/realloc rather than
Calloc/Realloc -- does it make sense to go through the code
replacing these all and see if that fixes the problem?

 cheers
   Ben Bolker

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] weights in lmer

[Chris O'Brien]

Dear R users,

Prof. Ripley just corrected my understanding of the use of weights in glm, 
which I thought would allow me to correctly use lmer.  However I'm still 
having problems.

My data takes the form of # of infected and uninfected individuals that 
were measured over time under different treatments.  I'm using lmer to 
adjust for the repeated measures over time.

In fitting the model:

 > model1=lmer(y~treatment+(time | ID),family=binomial,weights = total)

where y = proportion of animals infected (number infected/total)
 total = number of infected + number uninfected

this returns the error:

Error in lmer(y ~ treatment + (time | ID), family = binomial, weights = 
total) :
 object `weights' of incorrect type

Can anyone tell me what's wrong with the weights?

cheers,
Chris





Chris O'Brien
Sonoran Desert Research Station and
School of Natural Resources
Biological Sciences, room 125
University of Arizona
Tucson, AZ 85721
work (520) 623-3720
[EMAIL PROTECTED]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] labelling plots with ancillary data in data.frame

[Julian Burgos]

Hi Wesley,

You can use the text() function to add text to an existing plot.   See  
?text.

Julian

Wesley Roberts wrote:
> Hi All,
>
> I am busy using R to do some regression modelling and have been using 
> plot(x,y,"") to visualise my variables. I would now like to label my points 
> using data stored in the data.frame used for the regression analysis. For 
> example each of my data points is made up of a field measured forest volume 
> value and a remotely sensed vegetation estimate (NDVI). Each point is an 
> enumeration plot and I would like to label each the points in the 
> xy-scatterplot with their respective plot numbers. Is this possible in R, if 
> so how do I go about doing it?
>
> Many thanks for your help
>
> Wesley
>
> Wesley Roberts MSc.
> Researcher: Forest Assessment (Remote Sensing & GIS)
> Forestry and Forest Products Research Centre
> CSIR
> Tel: +27 (31) 242-2353
> Fax: +27 (31) 261-1216
> http://ffp.csir.co.za/
>
> "To know the road ahead, ask those coming back."
> - Chinese proverb
>
>
>   

-- 
Julian M. Burgos

Fisheries Acoustics Research Lab
School of Aquatic and Fishery Science
University of Washington

1122 NE Boat Street
Seattle, WA  98105 

Phone: 206-221-6864

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] weights in GAMs (package mgcv)

[Prof Brian Ripley]

Let's simplify to a linear model.  If your covariates have uncertainties, 
most likely a linear regression is not appropriate.  This sounds like an 
'errors in measurements' model, as covered in

@Book{Fuller.87,
   author   = "Fuller, Wayne A.",
   title= "Measurement Error Models",
   publisher= "John Wiley and Sons",
   address =  "New York",
   year = "1987",
   ISBN = "0-471-86187-1",
}

in which there is a true covariate that enters the model, but it is only 
observed with measurement error (or similar scenarios).

This is hard enough for linear models, without thinking about non-normal 
models or extensions beyond linear predictors.  The GLM (including GAM) 
estimation process assumes various things, including that the covariates 
that enter into the model are fixed (possibly by conditioning on them) and 
known.

On Tue, 14 Aug 2007, Julian Burgos wrote:

> Dear list,
>
> I?m using the ?mgcv? package to fit some GAMs. Some of my covariates are
> derived quantities and have an associated standard error, so I would
> like to incorporate this uncertainty into the GAM estimation process.
> Ideally, during the estimation process less importance would be given to
> observations whose covariates have high standard errors.
>
> The gam() function in the ?mgcv? package has a ?weights? argument.
> According to the package documentation, this can be used to provide
> prior weights to the data. This argument (as far as I understand) takes
> a vector of the same length of the data with numeric values higher than
> zero. So it seems that I should combine the standard errors of all
> covariates into a single vector and use it as weights. But it is not
> obvious to me how to do this, given that the covariates have different
> units and ranges of values.

Actually this is just taken from glm(), and case weights are part of the 
definition of a GLM.  In so far as I understand your scenario, you do not 
have a GLM.

> Is there any way to provide weights to the covariates directly (for
> example providing a matrix of n x m values, where n=number of covariates
> and m=number of observations)?
>
> Thanks,
>
> Julian
>
> Julian M. Burgos
>
> Fisheries Acoustics Research Lab
> School of Aquatic and Fishery Science
> University of Washington
>
> 1122 NE Boat Street
> Seattle, WA  98105
>
> Phone: 206-221-6864
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Import of Access data via RODBC changes column name ("NO" to "Expr1014") and the content of the column

[Prof Brian Ripley]

On Tue, 14 Aug 2007, Maciej Hoffman-Wecker wrote:


Dear Professor Ripley,

Thank you very much for your response. I send the problem, as I didn't 
have any more ideas were to search for the reason. I didn't say this is 
a R bug, knowing the responses on such mails.-)


But I succeeded in developing a tiny example, that reproduces the bug 
(wherever it is).


Thank you, that was helpful: much easier to follow that the previous code.

...


library(RODBC)
.con <- odbcConnectAccess("./test2.mdb")
(.d <- try(sqlQuery(.con, "select * from Tab1")))

 F1 NO F2
1  1  1  1
2  2  2  2
3  0 NA  1
4  1  0  0

(.d <- try(sqlQuery(.con, "select F1 , NO , F2 from Tab1")))

 F1 Expr1001 F2
1  10  1
2  20  2
3  00  1
4  10  0

close(.con)


So the problem occurs if the column names are specified within the query.
Is the query "select F1 , NO , F2 from Tab1" invalid?


I believe so. 'NO' is an SQL92 and ODBC reserved word, at least according 
to http://www.bairdgroup.com/reservedwords.cfm


See also http://support.microsoft.com/default.aspx?scid=kb;en-us;286335
which says

  For existing objects with names that contain reserved words, you can
  avoid errors by surrounding the object name with brackets ([ ]).

and lists 'NO' as a reserved word.  RODBC quotes all column names it uses 
to be sure (and knows about most non-standard quoting mechanisms from the 
ODBC driver in use).  But this was a query you generated and so you need 
to do the quoting.


Regarding the memory issue, I _knew_ that there must be a reason for the 
running out of memory space. Sorry for not being more specific. My 
question than is:


Is there a way to 'reset' the environment without quitting R and 
restarting it?


Sorry, no.  You cannot move objects in memory.

But why '477Mb' is coming up is still unexplained, and suggests that the 
machine has a peculiar amount of memory or some flag has been used.





Thank you for your help.

Kind regards,
Maciej


-Ursprüngliche Nachricht-
Von: Prof Brian Ripley [mailto:[EMAIL PROTECTED]
Gesendet: Dienstag, 14. August 2007 11:51
An: Maciej Hoffman-Wecker
Cc: r-help@stat.math.ethz.ch
Betreff: Re: [R] Import of Access data via RODBC changes column name ("NO" to 
"Expr1014") and the content of the column

On Tue, 14 Aug 2007, Maciej Hoffman-Wecker wrote:



Dear all,

I have some problems with importing data from an Access data base via
RODBC to R. The data base contains several tables, which all are
imported consecutively. One table has a column with column name "NO".
If I run the code attached on the bottom of the mail I get no
complain, but the column name (name of the respective vector of the
data.frame) is "Expr1014" instead of "NO". Additionally the original
column (type
"text") containes "0"s and missings, but the imported column contains
"0"s only (type "int"). If I change the column name in the Access data
base to "NOx", the import works fine with the right name and the same
data.

Previously I generated a tiny Access data base which reproduced the
problem. To be on the safe site I installed the latest version (2.5.1)
and now the example works fine, but within my production process the
error still remaines. An import into excel via ODBC works fine.

So there is no way to figure it out whether this is a bug or a
feature.-)


It's most likely an ODBC issue, but you have not provided a reproducible 
example.


The second problem I have is that when I rerun "rm(list = ls(all =
T)); gc()" and the import several times I get the following error:

Error in odbcTables(channel) : Calloc could not allocate (263168 of 1)
memory In addition: Warning messages:
1: Reached total allocation of 447Mb: see help(memory.size) in:
odbcQuery(channel, query, rows_at_time)
2: Reached total allocation of 447Mb: see help(memory.size) in:
odbcQuery(channel, query, rows_at_time)
3: Reached total allocation of 447Mb: see help(memory.size) in:
odbcTables(channel)
4: Reached total allocation of 447Mb: see help(memory.size) in:
odbcTables(channel)

which is surprising to me, as the first two statements should delete
all


How do you _know _what they 'should' do?  That only deletes all objects in the 
workspace, not all objects in R, and not all memory blocks used by R.

Please do read ?"Memory-limits" for the possible reasons.

Where did '447Mb' come from?  If this machine has less than 2Gb of RAM, buy 
some more.



objects and recover the memory. Is this only a matter of memory? Is
there any logging that reduces the memory? Or is this issue connected to
the upper problem?

I added the code on the bottom - maybe there is some kind of misuse I
lost sight of. Any hints are appreciated.

Kind regards,
Maciej


version

  _
platform   i386-pc-mingw32
arch   i386
os mingw32
system i386, mingw32
status
major  2
minor  5.1
year   2007
month  06
day27
svn rev42083
language   R
version.string R version 2

Re: [R] glm(family=binomial) and lmer

[Prof Brian Ripley]

On Tue, 14 Aug 2007, Chris O'Brien wrote:

> Dear R users,
>
> I've notice that there are two ways to conduct a binomial GLM with binomial
> counts using R.  The first way is outlined by Michael Crawley in his
> "Statistical Computing book" (p 520-521):

and in the places he got it from (it is not his original work).

These are not the only two ways, and they are not the same analyses as the 
saturated models differ.  The usual way to use weights is

y <- dead/batch
model3 <- glm(y ~ log(dose), binomial, weights=batch)
summary(model3)

and internally glm converts models with a two-column response to this 
form, for it is in this form the binomial fits into the GLM framework.

See the White Book or MASS (even the 1994 edition).


> >dose=c(1,3,10,30,100)
> >dead = c(2,10,40,96,98)
> >batch=c(100,90,98,100,100)
> >response = cbind(dead,batch-dead)
> >model1=glm(y~log(dose),binomial)
> >summary(model1)
>
> Which returns (in part):
> Coefficients:
> Estimate Std. Error z value Pr(>|z|)
> (Intercept)  -4.5318 0.4381  -10.35   <2e-16 ***
> log(dose) 1.9644 0.1750   11.22   <2e-16 ***
> Null deviance: 408.353  on 4  degrees of freedom
> Residual deviance:  10.828  on 3  degrees of freedom
> AIC: 32.287
>
> Another way to do the same analysis is to reformulate the data, and use GLM
> with weights:
>
> >y1=c(rep(0,5),rep(1,5))
> >dose1=rep(dose,2)
> >number = c(batch-dead,dead)
> >data1=as.data.frame(cbind (y1,dose,number))
> >model2=glm(y1~log(dose1),binomial,weights=number,data=data1)
> >summary(model2)
>
> Which returns:
>
> Coefficients:
> Estimate Std. Error z value Pr(>|z|)
> (Intercept)  -4.5318 0.4381  -10.35   <2e-16 ***
> log(dose1)1.9644 0.1750   11.22   <2e-16 ***
> (Dispersion parameter for binomial family taken to be 1)
> Null deviance: 676.48  on 9  degrees of freedom
> Residual deviance: 278.95  on 8  degrees of freedom
> AIC: 282.95
>
> Number of Fisher Scoring iterations: 6
>
> These two methods are similar in the parameter estimates and standard
> errors, however the deviances, their d.f., and AIC differ.  I take the
> first method to be the correct one.

This form has ten obeservations of groups with weights 2,98,10,80 

> However, I'm really interested in conducting a GLM binomial mixed model,
> and I am unable to figure out how to use the first method with the lmer
> function from the lme4 library, e.g.
>
> >model3=lmer(y~log(dose)+time|ID)# the above example data doesn't have
> the random effect, but my own data set does.
>
>   Does anyone have any suggestions?
>
> thanks,
> chris
>
> Thanks,
> Chris O'Brien
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] weights in GAMs (package mgcv)

[Julian Burgos]

Dear list,

I’m using the ‘mgcv’ package to fit some GAMs. Some of my covariates are 
derived quantities and have an associated standard error, so I would 
like to incorporate this uncertainty into the GAM estimation process. 
Ideally, during the estimation process less importance would be given to 
observations whose covariates have high standard errors.

The gam() function in the ‘mgcv’ package has a ‘weights’ argument. 
According to the package documentation, this can be used to provide 
prior weights to the data. This argument (as far as I understand) takes 
a vector of the same length of the data with numeric values higher than 
zero. So it seems that I should combine the standard errors of all 
covariates into a single vector and use it as weights. But it is not 
obvious to me how to do this, given that the covariates have different 
units and ranges of values.

Is there any way to provide weights to the covariates directly (for 
example providing a matrix of n x m values, where n=number of covariates 
and m=number of observations)?

Thanks,

Julian

Julian M. Burgos

Fisheries Acoustics Research Lab
School of Aquatic and Fishery Science
University of Washington

1122 NE Boat Street
Seattle, WA  98105 

Phone: 206-221-6864

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Import of Access data via RODBC changes column name ("NO" to "Expr1014") and the content of t he column

[Dieter Menne]

Maciej Hoffman-Wecker  bioskin.de> writes:

...
> But I succeeded in developing a tiny example, that reproduces the bug
(wherever it is).
> 
> I generated a small Access data base "test2.mdb" with one table "Tab1" and
following columns:
. 

> > library(RODBC)
> > .con <- odbcConnectAccess("./test2.mdb")
> > (.d <- try(sqlQuery(.con, "select * from Tab1")))

> > (.d <- try(sqlQuery(.con, "select F1 , NO , F2 from Tab1")))
>   F1 Expr1001 F2
> 1  10  1
> 2  20  2
> 3  00  1
> 4  10  0
> > close(.con)

NO is a reserved word in ODBC (or where...). Whenever you see Exprxxx in
columns, put the column name in []
d <- try(sqlQuery(con, "select F1 ,[NO], F2 from Tab1"))

works for me.

Dieter

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] labelling plots with ancillary data in data.frame

[Paul Hiemstra]

Hi Wesley,

Try the text() function. An example:

a = rep(10,10)
b = seq(1,10)
plot(a,b)
text(a,b, labels = b, pos = 4, offset = 0.7)
?text

hth,

Paul

Wesley Roberts schreef:
> Hi All,
>
> I am busy using R to do some regression modelling and have been using 
> plot(x,y,"") to visualise my variables. I would now like to label my points 
> using data stored in the data.frame used for the regression analysis. For 
> example each of my data points is made up of a field measured forest volume 
> value and a remotely sensed vegetation estimate (NDVI). Each point is an 
> enumeration plot and I would like to label each the points in the 
> xy-scatterplot with their respective plot numbers. Is this possible in R, if 
> so how do I go about doing it?
>
> Many thanks for your help
>
> Wesley
>
> Wesley Roberts MSc.
> Researcher: Forest Assessment (Remote Sensing & GIS)
> Forestry and Forest Products Research Centre
> CSIR
> Tel: +27 (31) 242-2353
> Fax: +27 (31) 261-1216
> http://ffp.csir.co.za/
>
> "To know the road ahead, ask those coming back."
> - Chinese proverb
>
>
>   


-- 
Drs. Paul Hiemstra
Department of Physical Geography
Faculty of Geosciences
University of Utrecht
Heidelberglaan 2
P.O. Box 80.115
3508 TC Utrecht
Phone:  +31302535773
Fax:+31302531145
http://intamap.geo.uu.nl/~paul

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] glm(family=binomial) and lmer

[Chris O'Brien]

Dear R users,

I've notice that there are two ways to conduct a binomial GLM with binomial 
counts using R.  The first way is outlined by Michael Crawley in his 
"Statistical Computing book" (p 520-521):

 >dose=c(1,3,10,30,100)
 >dead = c(2,10,40,96,98)
 >batch=c(100,90,98,100,100)
 >response = cbind(dead,batch-dead)
 >model1=glm(y~log(dose),binomial)
 >summary(model1)

Which returns (in part):
Coefficients:
 Estimate Std. Error z value Pr(>|z|)
(Intercept)  -4.5318 0.4381  -10.35   <2e-16 ***
log(dose) 1.9644 0.1750   11.22   <2e-16 ***
Null deviance: 408.353  on 4  degrees of freedom
Residual deviance:  10.828  on 3  degrees of freedom
AIC: 32.287

Another way to do the same analysis is to reformulate the data, and use GLM 
with weights:

 >y1=c(rep(0,5),rep(1,5))
 >dose1=rep(dose,2)
 >number = c(batch-dead,dead)
 >data1=as.data.frame(cbind (y1,dose,number))
 >model2=glm(y1~log(dose1),binomial,weights=number,data=data1)
 >summary(model2)

Which returns:

Coefficients:
 Estimate Std. Error z value Pr(>|z|)
(Intercept)  -4.5318 0.4381  -10.35   <2e-16 ***
log(dose1)1.9644 0.1750   11.22   <2e-16 ***
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 676.48  on 9  degrees of freedom
Residual deviance: 278.95  on 8  degrees of freedom
AIC: 282.95

Number of Fisher Scoring iterations: 6

These two methods are similar in the parameter estimates and standard 
errors, however the deviances, their d.f., and AIC differ.  I take the 
first method to be the correct one.
However, I'm really interested in conducting a GLM binomial mixed model, 
and I am unable to figure out how to use the first method with the lmer 
function from the lme4 library, e.g.

 >model3=lmer(y~log(dose)+time|ID)# the above example data doesn't have 
the random effect, but my own data set does.

   Does anyone have any suggestions?

thanks,
chris

Thanks,
Chris O'Brien

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Problem with "by": does not work with ttest (but with lme)

[Greg Snow]

I think at least part of your problem is that in the lm example you use
data=x (correct), but in the t.test example you use data=warpbreaks, so
in that case it is uing the full data set, not just the portion passed
by the 'by' function.  Try the t.test example again with data=x and see
what happens.

Hope this helps,

-- 
Gregory (Greg) L. Snow Ph.D.
Statistical Data Center
Intermountain Healthcare
[EMAIL PROTECTED]
(801) 408-8111
 
 

> -Original Message-
> From: [EMAIL PROTECTED] 
> [mailto:[EMAIL PROTECTED] On Behalf Of Daniel Stahl
> Sent: Tuesday, August 14, 2007 7:11 AM
> To: r-help@stat.math.ethz.ch
> Subject: [R] Problem with "by": does not work with ttest (but 
> with lme)
> 
> Hello,
> 
> I would like to do a large number of e.g. 1000 paired ttest 
> using the by-function. But instead of using only the data 
> within the 1000 groups, R caclulates 1000 times the ttest for 
> the full data set(The same happens with Wilcoxon test). 
> However, the by-function works fine with the lme function.
> Did I just miss something or is it really not working? If 
> not, is there any other possibility to avoid loops? 
> Thanks
> Daniel
> 
> Here is the R help example for "by" 
>  require(stats)
>  attach(warpbreaks)
>  by(warpbreaks, tension, function(x) lm(breaks ~ wool, data = 
> x)) *->works great 
> by(warpbreaks,tension,function(x)t.test(breaks ~ 
> wool,data=warpbreaks,paired = TRUE)) *Same output for each 
> level of tension:
> 
> tension: L
> 
>   Paired t-test
> 
> data:  breaks by wool
> t = 1.9956, df = 26, p-value = 0.05656
> alternative hypothesis: true difference in means is not equal to 0
> 95 percent confidence interval:
> -0.1735803 11.7291358
> sample estimates:
> mean of the differences
> 5.78
> 
> --
> --
> 
> tension: M
> 
>   Paired t-test
> 
> data:  breaks by wool
> t = 1.9956, df = 26, p-value = 0.05656
> alternative hypothesis: true difference in means is not equal to 0
> 95 percent confidence interval:
> -0.1735803 11.7291358
> sample estimates:
> mean of the differences
> 5.78
> 
> --
> --
> 
> tension: H
> 
>   Paired t-test
> 
> data:  breaks by wool
> t = 1.9956, df = 26, p-value = 0.05656
> alternative hypothesis: true difference in means is not equal to 0
> 95 percent confidence interval:
> -0.1735803 11.7291358
> sample estimates:
> mean of the differences
> 5.78
> 
> 
> 
> 
> 
> 
> 
> --
> 
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide 
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] kernlab ksvm() cross-validation prediction response vector

[strinz]

Hello,

I would like to know, whether for the support vector classification function 
ksvm()
the response values stored in [EMAIL PROTECTED] are cross validated 
outputs/predictions:

Example code from package kernlab, function ksvm:
library(kernlab)
## train a support vector machine
filter <- 
ksvm(type~.,data=spam,kernel="rbfdot",kpar=list(sigma=0.05),C=5,cross=3)
filter
[EMAIL PROTECTED]

if not:
what is the easiest way to obtain x-fold cross validated
predicted values of the instances of a data set ?
Does it have to be implemented by oneself?

thanks a lot again
best regards 
Björn

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Limma - 2x2 factorial design matrix

[tee` suesean]

Hi all,
 
I'm working on microarray and currently analyzing the microrarray data 
using limmaGUI. Loop design has been applied in this experiment. This is a 2X2 
factorial experiment where there are control and treatment at 2 different time 
points, week 6 and 9. The experimental design is almost the same as the 
limmaGUI work example: Weaver Data set. I would like to look at the effect of 
time, treatment and interaction. I have tried to follow the limmaGUI help 
example to construct the design matrix. The target file and design matrix are 
as follow: (C - control, M - treated)
 
SlideNumberFileNameCy3Cy5
M6.C6tss37-normnew.gprM6C6
C6.C9tss5-new.gprC6C9
C9.M9tss9-new.gprC9M9
M9.M6tss24-normnew.gprM9M6
C6.M6tss31-normnew.gprC6M6
M6.M9tss2-normnew.gprM6M9
M9.C9tss12-normnew.gprM9C9
C9.C6tss25-normnew.gprC9C6

 
SlideNumberTreatmentTimeInteraction
M6.C6-100
C6.C9010
C9.M9101
M9.M60-1-1
C6.M6100
M6.M9011
M9.C9-10-1
C9.C60-10

 
According to the design matrix above, may I know what is the meaning of time 
effect? Does it look at the differential expression of Control sample at 
different time? How about treatment effect? Is it refer to the differential 
gene expression for Treated sample at week 6? For the effect of interaction, is 
it refer to both time and treatment refer, which is refer to M9?
 
Besides that, I can not figure out the second example of design matrix, Mt for 
Weaver data set in LimmaGUI help file. Is the design matrix looking at the 
effect of different time point for mutant sample? Mt11 and Mt21 refer to the 
effect of mutant at day 11 and day 21 alone? What is the logically theory 
behind to asign 1 or -1 to slides P21WT.P11WT for effect of Mt11 and Mt21? Does 
it has anything to do with the vector? Or we will subtract it with Wt at both 
time point in order for us to look at the effect of Mt11 and Mt21?
 
I use Lucidea Universal ScoreCard as my spike in control. Do I need to flag 
them out before I carry out the analysis? Although I have include a Spot Type 
file, LimmaGUI still read all the spots as genes. 
 
Thank you 
 
 
Warm regards,
Tee Sue Sean
Department of Cell and Molecular Biology,
Faculty of Biotechnology and Biomolecular Sciences,
University Putra Malaysia, 
43400 Serdang, Selangor, Malaysia.



Choose the right car based on your needs. Check out Yahoo! Autos new Car Finder 
tool.


   




[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Using sunflowerplot to add points in a xyplot panel

[Dieter Menne]

Ronaldo Reis Junior  gmail.com> writes:

> 
> Hi,
> 
> I use panel.points to add points to a xyplot graphic. But I like to use the 
> sunflowerplot to plot my points because this is very superimposed. It is 
> possible to use this? I try but it dont work directly. It may be need to put 
> this function inside a panel.??? 

sunflower is "old style R graphics", and xyplot is lattice/trellis. these two
normally don't get along together too well. Paul Murrell's grid which is the
base of lattice, can do some magic to bring them together.

Dieter


## Modfied after Cris Bergstresser and Paul Murrell
##http://finzi.psych.upenn.edu/R/Rhelp02a/archive/73780.html
library(grid);
library(lattice);

sunpanel <- function(x, y, subscripts, ...) {
  pushViewport(viewport(x = 0.5, y = 0.5, just = "center"));
  sunflowerplot(x, y, axes = FALSE, xlab = "", ylab = "");
  popViewport(1);
} 

x = round(runif(100, 1, 5));
y = round(runif(100, 1, 5));
print(xyplot(y ~ x, aspect = "iso", panel = sunpanel));

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Problem with

[Dieter Menne]

Daniel Stahl  operamail.com> writes:

> 
> I would like to do a large number of e.g. 1000 paired ttest using the
by-function.

Technical stuff aside: do you really want to do 1000 paired t-tests?
Followed by a *** pick-nic?

Mmm...

Dieter

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] invert 160000x160000 matrix

[Douglas Bates]

On 8/14/07, Jiao Yang <[EMAIL PROTECTED]> wrote:

> Tir,

> Thank you very much for the note.  I'm using an algorithm to analyze a data 
> set of 400 variables and the algorithm uses the inverse of a 16x16 
> positive definite matrix.  The matrix is stored as a text file.

A meta-theorem in numerical linear algebra is that if your algorithm
involves finding the inverse of a matrix then you need a better
algorithm.  When the matrix is of size p^2 by p^2 for p variables and
you say that you need to invert it you really need a better algorithm.

You need to do a lot of reading on numerical linear algebra and a lot
of thinking about the algorithm before you can proceed.  The good news
is that you will have a lot of time to do this while you are waiting
from the result from any method that you could possibly try to do the
calculation the way you have phrased it.


> Can you please give some reference to  in-place inverse technique?  I googled 
> it, yet did not get any good results.  Thank you for your time.
>
> Best regards,
> Jiao
>
>
> - Original Message -
> From: "Patnaik, Tirthankar " <[EMAIL PROTECTED]>
> Date: Tuesday, August 14, 2007 2:07 am
> Subject: RE: [R] invert 16x16 matrix
> To: Moshe Olshansky <[EMAIL PROTECTED]>, Paul Gilbert <[EMAIL PROTECTED]>, 
> Jiao Yang <[EMAIL PROTECTED]>
> Cc: r-help@stat.math.ethz.ch
>
>
> > A variety of tricks would need to be used to invert a matrix of this
> > size. If there are any other properties of the matrix that you know
> > (symmetric, positive definite, etc, sparse) then they could be useful
> > too. You could partition the matrix first, then use an in-place
> > inverse technique for each block to individually calculate the
> > blocks-inverses, then combine to get the inverse of the initial
> > matrix. Again, if the implementation is actually solving an Ax-B = 0
> > system of equations, then there are specific methods for these too,
> > like an LU decomp, for instance. You might also want to check out some
> > texts for this, like the Numerical Recipes.
> > How's the matrix stored right now?
> >
> > Best,
> > -Tir
> >
> >
> >
> > Tirthankar Patnaik
> > India Strategy
> > Citigroup Investment Research
> > +91-22-6631 9887
> >
> > For important disclosures regarding Citigroup Investment Research,
> > including with respect to any issuers mentioned herein, please refer
> > to the Citigroup Investment Research disclosure website at
> >
> >
> > > -Original Message-
> > > From: [EMAIL PROTECTED]
> > > [ On Behalf Of Moshe Olshansky
> > > Sent: Tuesday, August 14, 2007 6:40 AM
> > > To: Paul Gilbert; Jiao Yang
> > > Cc: r-help@stat.math.ethz.ch
> > > Subject: Re: [R] invert 16x16 matrix
> > >
> > > While inverting the matrix may be a problem, if you need to
> > > solve an equation A*x = b you do not need to invert A, there
> > > exist iterative methods which do need A or inv(A) - all you
> > > need to provide is a function that computes A*x for an
> > > arbitrary vector x.
> > > For such a large matrix this may be slow but possible.
> > >
> > > --- Paul Gilbert <[EMAIL PROTECTED]>
> > > wrote:
> > >
> > > > I don't think you can define a matrix this large in R, even if you
> >
> > > > have the memory. Then, of course, inverting it there may be other
> >
> > > > programs that have limitations.
> > > >
> > > > Paul
> > > >
> > > > Jiao Yang wrote:
> > > > > Can R invert a 16x16 matrix with all
> > > > positive numbers?  Thanks a lot!
> > > > >
> > > > > __
> > > > > R-help@stat.math.ethz.ch mailing list
> > > > >
> > > > > PLEASE do read the posting guide
> > > >
> > > > > and provide commented, minimal, self-contained,
> > > > reproducible code.
> > > >
> > > ==
> > > ==
> > > >
> > > > La version française suit le texte anglais.
> > > >
> > > >
> > > --
> > > --
> > > >
> > > > This email may contain privileged and/or confidential
> > > > inform...{{dropped}}
> > > >
> > > > __
> > > > R-help@stat.math.ethz.ch mailing list
> > > >
> > > > PLEASE do read the posting guide
> > > >
> > > > and provide commented, minimal, self-contained, reproducible code.
> > > >
> > >
> > > __
> > > R-help@stat.math.ethz.ch mailing list
> > >
> > > PLEASE do read the posting guide
> > >
> > > and provide commented, minimal, self-contained, reproducible code.
> > >
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.c

Re: [R] diffusing GIS data in maps

[Paul Hiemstra]

Hi Lawrence,

You could use the gstat (geostatistics) package to perform an 
interpolation, it also requires the package sp. It offers inverse 
distance interpolation and several forms of kriging. Making an 
interpolation would look something like:

library(gstat) # Also loads sp
coordinates(geocode) = ~LAT + LONG  # make a spatial object out of 
geocode
grid = spsample(geocode, type = "regular", n = 1000)   # make a map of 
target locations.
# Regular grid. n is the number 
of cells, increase this for more detail in the resulting map
gridded(grid) = TRUE # Identifying this object as a grid
interpolated = krige(VALUE ~ 1, geocode, grid)   # Do inverse distance 
interpolation, data is 'geocode', target locations are 'grid'.
spplot(interpolated, "var1.pred", sp.layout = list("sp.points", 
geocode))   # Plot the predictions (grid) and the points (geocode)

hope this helps,

Paul

Lawrence D. Brenninkmeyer schreef:
> Hi-
>
> I am trying to find a way to diffuse GIS data on a European map. I have a
> dataset consisting of particular locations scattered across Europe,
> along with magnitude and value information. I can plot these as discrete
> points with something like the following:
>
> "geocode" is a dataframe with four columns: LAT; LONG; MAGNITUDE;VALUE.
>
> library(maps)
> library(mapdata)
> map("worldHires", regions=c("Germany", "Belgium", "Netherlands"))
> points(geocode$LONG, geocode$LAT, cex=geocode$MAGNITUDE / 2500,
> col=rainbow(length(geocode$VALUE), start=0, end=.4)[rank(geocode$VALUE)])
>
> This gives me a map of Europe with my datapoints highlighted in two ways:
> magnitude is represented by the size of the point, and value is
> represented by the color.
>
> However, what I would really like is for there to be some sort of
> diffusion, such that instead of discrete points, the European map is
> covered in color so I can see more clearly whether there are regional
> patterns (something that will presumably look like this contour chart:
>   only
> on the European map).
>
> I have absolutely no idea where to start because I can't find a function
> that will allow me to diffuse the datapoints on a map.
>
> thank you for any help
> ldb
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>   


-- 
Drs. Paul Hiemstra
Department of Physical Geography
Faculty of Geosciences
University of Utrecht
Heidelberglaan 2
P.O. Box 80.115
3508 TC Utrecht
Phone:  +31302535773
Fax:+31302531145
http://intamap.geo.uu.nl/~paul

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] {grid} plain units with non NULL data arguments

[Prof Brian Ripley]

On Tue, 14 Aug 2007, Wolfram Fischer wrote:

> In help(unit) I read:
>
> The 'data' argument must be a list when the 'unit.length()'
> is greater than 1.  For example, 'unit(rep(1, 3), c("npc",
> "strwidth", "inches"), data=list(NULL, "my string", NULL))'.
>
> In the newest R-versions it is not anymore allowed to let strings
> in the data-argument for plain units, otherwise one gets the
> following error:
>Non-NULL value supplied for plain unit
>
> I have some labels. Between them I wanted to set a distance of 1.5 lines.
> (I wanted to use that for a grid.layout for a legend:
> The space is for the symbols.)
>
>labels <- c( ':', 'a', 'bb', 'ccc', '', 'e' )
>n <- length( labels )
>s <- as.list( c( labels[1], rep( labels[-1], each=2 ) ) )
>u <- unit( data=s, x=c( 1, rep( c( 1.5, 1 ), n-1 ) ),
>units=c( 'strwidth', rep( c( 'lines', 'strwidth' ), n-1 ) ) )
>
> How can I insert the NULL values into the list ``s''?
>
> To fill every second element of s with NULL, I tried:
>s[ 2 * ( 1 : length( labels[-1] ) ) ] <- NULL
> But this deletes every second element.

A value of list(NULL) is correct for inserting NULLs into lists.
(More generally to substitute in a list you need a list value.)

> The following would work:
>s[ 2 * ( 1 : length( labels[-1] ) ) ] <- NA
> But unit() does not accept NAs.

More to the point, it does not accept logical vectors as NULL values.

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Question about unicode characters in tcltk

[Peter Dalgaard]

R Help wrote:
> hello list,
>
> Can someone help me figure out why the following code doesn't work?
> I'm trying to but both Greek letters and subscripts into a tcltk menu.
>  The code creates all the mu's, and the 1 and 2 subscripts, but it
> won't create the 0.  Is there a certain set of characters that R won't
> recognize the unicode for?  Or am I input the \u2080 incorrectly?
>
> library(tcltk)
> m <-tktoplevel()
> frame1 <- tkframe(m)
> frame2 <- tkframe(m)
> frame3 <- tkframe(m)
> entry1 <- tkentry(frame1,width=5,bg='white')
> entry2 <- tkentry(frame2,width=5,bg='white')
> entry3 <- tkentry(frame3,width=5,bg='white')
>
> tkpack(tklabel(frame1,text='\u03bc\u2080'),side='left')
> tkpack(tklabel(frame2,text='\u03bc\u2081'),side='left')
> tkpack(tklabel(frame3,text='\u03bc\u2082'),side='left')
>
> tkpack(frame1,entry1,side='top')
> tkpack(frame2,entry2,side='top')
> tkpack(frame3,entry3,side='top')
>
>   
Odd, but I think not an R issue. I get weirdness in wish too. Try this

% toplevel .a
.a
% label .a.b -text \u03bc\u2080 -font {Roman -10}
.a.b
% pack .a.b
% .a.b configure
{-activebackground
[]
{-text text Text {} μ₀} {-textvariable textVariable Variable {} {}}
{-underline underline Underline -1 -1} {-width width Width 0 0}
{-wraplength wrapLength WrapLength 0 0}
% .a.b configure -font {Helvetica -12 bold} # the default, now shows \u2080
% .a.b configure -font {Roman -10} # back to Roman, *still* shows \u2080

???!!!

-- 
   O__   Peter Dalgaard Øster Farimagsgade 5, Entr.B
  c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K
 (*) \(*) -- University of Copenhagen   Denmark  Ph:  (+45) 35327918
~~ - ([EMAIL PROTECTED])  FAX: (+45) 35327907

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Linear Regression with slope equals 0

[ONKELINX, Thierry]

Dear Ed,

In my opinion you don't need to set the slope to zero. Just test if the
slope in lm(d ~ t) is significant. If it is significant then you have
evidence that the slope is NOT zero. But when it is not significant (and
in your example it is), you can't say that it is zero. But doing some
power calculations allows you to estimate the smallest detectable slope.
The estimated slope in your example is 0.49. Your design has a power of
13% to detect slopes of this size.  The real slope has to be about 1.75
in order to have 80% power. Let's say that it would matter if the slope
is about 2 or larger and it won't matter is it is below 2. The power to
detect a slope of 2 is 89%, hence you would likely get a significant
slope if it was larger then 2. Since the slope was not significant, it's
safe to say that the slope is not larger then 2.
But if a slope of 0.5 would matter, you couldn't make this kind of
assumptions because the power to detect a slope of 0.5 is only 13%. So
the change of rejecting the null hypothesis (slope = 0) when slope = 0.5
is too small.

HTH,

Thierry

> power.trend <- function(repetitions = 5, x = c(0, 1), sd = 1, slope =
1, alpha = 0.05){
+   X <- rep(x, repetitions)
+   ncp <- slope ^ 2 * sum((X - mean(X))^2) / sd ^ 2
+   return(1 - pf(qf(1 - alpha, 1, length(X) - 2), 1, length(X) - 2, ncp
= ncp))
+ }
> 
> df <- data.frame(t = 1:6, d = c(303, 302, 304, 306, 307, 303))
> fit <- lm(d ~ t, data = df)
> summary(fit)

Call:
lm(formula = d ~ t, data = df)

Residuals:
   123456 
 0.04762 -1.43810  0.07619  1.59048  2.10476 -2.38095 

Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 302.4667 1.7849  169.45 7.28e-09 ***
t 0.4857 0.45831.060.349
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Residual standard error: 1.917 on 4 degrees of freedom
Multiple R-Squared: 0.2192, Adjusted R-squared: 0.02402 
F-statistic: 1.123 on 1 and 4 DF,  p-value: 0.349 

> power.trend(repetitions = 1, x = df$t, sd = sd(df$d), slope =
coef(fit)["t"])
[1] 0.1292668
> power.trend(repetitions = 1, x = df$t, sd = sd(df$d), slope = 1.75)
[1] 0.8021454




ir. Thierry Onkelinx
Instituut voor natuur- en bosonderzoek / Research Institute for Nature
and Forest
Cel biometrie, methodologie en kwaliteitszorg / Section biometrics,
methodology and quality assurance
Gaverstraat 4
9500 Geraardsbergen
Belgium
tel. + 32 54/436 185
[EMAIL PROTECTED]
www.inbo.be 

Do not put your faith in what statistics say until you have carefully
considered what they do not say.  ~William W. Watt
A statistical analysis, properly conducted, is a delicate dissection of
uncertainties, a surgery of suppositions. ~M.J.Moroney

 

> -Oorspronkelijk bericht-
> Van: [EMAIL PROTECTED] 
> [mailto:[EMAIL PROTECTED] Namens 
> [EMAIL PROTECTED]
> Verzonden: dinsdag 14 augustus 2007 13:37
> Aan: r-help@stat.math.ethz.ch
> Onderwerp: [R] Linear Regression with slope equals 0
> 
> 
>  Hi there, am trying to run a linear regression with a slope of 0.
> 
>  I have a dataset as follows
> 
>  t d
>  1 303
>  2 302
>  3 304
>  4 306
>  5 307
>  6 303
> 
>  I would like to test the significance that these points 
> would lie on a horizontal straight line.
> 
>  The standard regression lm(d~t) doesn't seem to allow the 
> slope to be set.
> 
>  Any help very welcome.
> 
>  ed
> 
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide 
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Comparing long species lists via Sorenson’s dissimilarity

[Richard Price]

I have 4 very large species lists and I would like to compare them.
   
  I have the following results from running Sorenson’s dissimilarity tests:
   
  Norfolk Fens compared to Suffolk Coastal Fens:
QS=0.583961142689298
  Norfolk Fens compared to Breckland Edge Fens:
QS=0.714896020281379
  Norfolk Fens compared to Other Fens: 
QS=0.78572348898302
   
  Suffolk Coastal Fens compared to Breckland Edge Fens:QS=0.78572348898302
  Suffolk Coastal Fens compared to Other Fens:  
QS=0.855011155835705
   
  Breckland Edge Fens compared to Other Fens:  
QS=0.175091076893185
   
  Does anyone know how to interpret the above results? I have read a paper that 
states that when comparing two samples if QS is less than 50% than the samples 
are considered dissimilar, what is 100%? I mean 50% of what?
   
  There are a number of other tests that I can run using R. These are squared 
euclidan, manhattan, bray-curtis, jaccard, ruzicka, (dis)similarity ratio, 
ochiai, cosine compliment, and Raup-crick. Would it be advantageous to run 
these now that I have my sorensons result? It is especially easy for me to run 
the Bray-Curtis all I would need to do is change the terms from ‘binary’ to 
‘minimum’ and re-run. Would Bray-Curtis be a Percentage Similarity unlike the 
others?
   
  Here is an example of what happens when I run sorensons:
   
  #Compare Breckland with Other
  > a=table(breckland)
  > J=sum(breckland*other)
  > A=sum(breckland^2)
  > B=sum(other^2)
   
  > brecklandOther <- designdist(a,method=(A+B-2*J)/(A+B), terms = c("binary"))
  [1] 0.1750911
  attr(,"call")
  designdist(x = a, method = (A + B - 2 * J)/(A + B), terms = c("binary"))
  attr(,"method")
  [1] "0.175091076893185 binary"
   
  I do not know why the following lines are included (perhaps R is trying to 
calculate the header row? But it doesn’t seem to stop the result getting out.
   
  attr(,"call")
  designdist(x = a, method = (A + B - 2 * J)/(A + B), terms = c("binary"))
  attr(,"method")
   
  I’ve also managed to get diversity indices from R but someone told me that 
comparing diversity indices is not good because of the way that they are 
calculated.
   
  >norfolkShannon
  3.653768 
  > suffolkShannon
   Suffolk 
  3.646138 
  > brecklandShannon
  Breckland 
   4.051742 
  > otherShannon
 Other 
  3.587403
   
  For the above statistical methods is there a way of laying out the data 
perhaps via a graph that will show how similar/dissimilar the sites are?
   
  Many thanks for yesterdays help and in advance for any help given today.
   
  Richard Price
  University of Birmingham.
   
   
   

[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Using sunflowerplot to add points in a xyplot panel

[Sundar Dorai-Raj]

Ronaldo Reis Junior said the following on 8/14/2007 7:08 AM:
> Hi,
> 
> I use panel.points to add points to a xyplot graphic. But I like to use the 
> sunflowerplot to plot my points because this is very superimposed. It is 
> possible to use this? I try but it dont work directly. It may be need to put 
> this function inside a panel.??? 
> 
> Thanks
> Ronaldo

You'll need to write your own panel function. Here's one shot at it. 
Most of the code is from ?sunflowerplot with added touches for lattice 
capability.

HTH,

--sundar

panel.sunflowerplot <- function(x, y, number, log = "", digits = 6, 
rotate = FALSE,
 cex.fact = 1.5, size = 1/8, seg.col = 
2, seg.lwd = 1.5, ...) {
   n <- length(x)
   if(missing(number)) {
 x <- signif(x, digits = digits)
 y <- signif(y, digits = digits)
 orderxy <- order(x, y)
 x <- x[orderxy]
 y <- y[orderxy]
 first <- c(TRUE, (x[-1] != x[-n]) | (y[-1] != y[-n]))
 x <- x[first]
 y <- y[first]
 number <- diff(c((1:n)[first], n + 1))
   } else {
 if(length(number) != n)
   stop("'number' must have same length as 'x' and 'y'")
 np <- number > 0
 x <- x[np]
 y <- y[np]
 number <- number[np]
   }
   n <- length(x)
   n.is1 <- number == 1
   cex <- trellis.par.get("plot.symbol")$cex
   if(any(n.is1))
 lpoints(x[n.is1], y[n.is1], cex = cex, ...)
   if(any(!n.is1)) {
 lpoints(x[!n.is1], y[!n.is1], cex = cex/cex.fact, ...)
 i.multi <- (1:n)[number > 1]
 ppin <- par("pin")
 pusr <- unlist(current.panel.limits())
 xr <- size * abs(pusr[2] - pusr[1])/ppin[1]
 yr <- size * abs(pusr[4] - pusr[3])/ppin[2]
 i.rep <- rep.int(i.multi, number[number > 1])
 z <- numeric()
 for (i in i.multi) z <- c(z, 1:number[i] + if (rotate) 
stats::runif(1) else 0)
 deg <- (2 * pi * z)/number[i.rep]
 lsegments(x[i.rep], y[i.rep], x[i.rep] + xr * sin(deg),
   y[i.rep] + yr * cos(deg), col = seg.col, lwd = seg.lwd)
   }
}

library(lattice)
xyplot(Petal.Width ~ Petal.Length, iris, panel = panel.sunflowerplot)

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] {grid} plain units with non NULL data arguments

[Wolfram Fischer]

In help(unit) I read:

 The 'data' argument must be a list when the 'unit.length()'
 is greater than 1.  For example, 'unit(rep(1, 3), c("npc",
 "strwidth", "inches"), data=list(NULL, "my string", NULL))'.

In the newest R-versions it is not anymore allowed to let strings
in the data-argument for plain units, otherwise one gets the
following error:
Non-NULL value supplied for plain unit

I have some labels. Between them I wanted to set a distance of 1.5 lines.
(I wanted to use that for a grid.layout for a legend:
The space is for the symbols.)

labels <- c( ':', 'a', 'bb', 'ccc', '', 'e' )
n <- length( labels )
s <- as.list( c( labels[1], rep( labels[-1], each=2 ) ) )
u <- unit( data=s, x=c( 1, rep( c( 1.5, 1 ), n-1 ) ),
units=c( 'strwidth', rep( c( 'lines', 'strwidth' ), n-1 ) ) )

How can I insert the NULL values into the list ``s''?

To fill every second element of s with NULL, I tried:
s[ 2 * ( 1 : length( labels[-1] ) ) ] <- NULL
But this deletes every second element.

The following would work:
s[ 2 * ( 1 : length( labels[-1] ) ) ] <- NA
But unit() does not accept NAs.


Regards - Wolfram

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] State Space Modelling

[Bernardo Ribeiro]

Here is one of the codes.

Thanks a lot...


module(finmetrics)

genvasicek.ssf = function(param, tau=NULL, freq=1/12){

## 1. Check for valid inputs

if(length(param) < 5)
stop("Parameters must have length greater than 4")
N = length(param) - 13
if (length(tau) != N)
stop("Length of Tau is inconsistent with Parameters")

## 2.Extract Parameters and Impose Constraints

k1 = (param[1])
k2 = (param[2])
k3 = (param[3])

lambda1 = param[4]
lambda2 = param[5]
lambda3 = param[6]

Delta = (param[7])

s1 = (param[8])
s2 = (param[9])
s3 = (param[10])

r21 = param[11]
r31 = param[12]
r32 = param[13]

st1 = (param[14])
st2 = (param[15])
st3 = (param[16])
st4 = (param[17])

## 3. Compute Matrixes

A1 = function(x,tau){
(1-exp(-x*tau))/x
}

u1 = -A1(k1,tau)  ## ok
u2 = -A1(k2,tau)  ## ok
u3 = -A1(k3,tau)  ## ok
H = cbind (-u1/tau,-u2/tau,-u3/tau)  ## ok

v1 = (lambda1/k1)*(tau - A1(k1,tau)) + (lambda2/k2)*(tau - A1(k2,tau)) +
(lambda3/k3)*(tau - A1(k3,tau))
v2 = +3 * Delta * tau  # Teste

v311 = ((s1*s1)/(k1*k1)) * (tau - A1(k1,tau) - A1(k1,tau) +
A1(2*k1,tau))
v312 = ((s1*s2*r21)/(k1*k2)) * (tau - A1(k1,tau) - A1(k2,tau) +
A1(k1+k2,tau))
v313 = ((s1*s3*r31)/(k1*k3)) * (tau - A1(k1,tau) - A1(k3,tau) +
A1(k1+k3,tau))

v322 = ((s2*s2)/(k2*k2)) * (tau - A1(k2,tau) - A1(k2,tau) +
A1(2*k2,tau))
v321 = ((s1*s2*r21)/(k1*k2)) * (tau - A1(k1,tau) - A1(k2,tau) +
A1(k1+k2,tau))
v323 = ((s2*s3*r32)/(k2*k3)) * (tau - A1(k2,tau) - A1(k3,tau) +
A1(k2+k3,tau))

v333 = ((s3*s3)/(k3*k3)) * (tau - A1(k3,tau) - A1(k3,tau) +
A1(2*k3,tau))
v331 = ((s1*s3*r31)/(k1*k3)) * (tau - A1(k3,tau) - A1(k1,tau) +
A1(k3+k1,tau))
v323 = ((s2*s3*r32)/(k2*k3)) * (tau - A1(k2,tau) - A1(k3,tau) +
A1(k2+k3,tau))

v = v1+v2+0.5*(v311+v312+v313+v322+v321+v323+v333+v331+v323)

d = -v

A = matrix(c(0,0,0,0,0,0,0,0,0),ncol=3)
A[1,1] = 1-k1*freq
A[1,2] = 0
A[1,3] = 0
A[2,1] = 0
A[2,2] = 1-k2*freq
A[2,3] = 0
A[3,1] = 0
A[3,2] = 0
A[3,3] = 1-k3*freq

ct = c(-lambda1*freq,-lambda2*freq,-lambda3*freq)

Q = matrix(c(0,0,0,0,0,0,0,0,0),ncol=3)
Q[1,1] = s1^2
Q[1,2] = s1*s2*r21
Q[1,3] = s1*s3*r31
Q[2,1] = s2*s1*r21
Q[2,2] = s2^2
Q[2,3] = s2*s3*r32
Q[3,1] = s1*s3*r31
Q[3,2] = s2*s3*r32
Q[3,3] = s3^2

## 4. Compute the State Space Form

mDelta = matrix(c(ct,d),ncol=1)

mPhi = rbind(A,H)

mOmega1 = matrix(c(Q,0,0,0,0,0,0,0,0,0,0,0,0),ncol=7)
mOmega2 =
matrix(c(0,0,0,0,0,0,0,0,0,0,0,0,diag(c(st1,st2,st3,st4))),ncol=7)

mOmega = rbind(mOmega1,mOmega2)

mSigma = matrix(c(0,0,0,0,0,0,0,0,0,0,0,0),ncol=3)

mSigma[1,1] = ((s1^2)/(2*k1))*(1-exp(-2*k1*freq))
mSigma[2,2] = ((s2^2)/(2*k2))*(1-exp(-2*k2*freq))
mSigma[3,3] = ((s3^2)/(2*k3))*(1-exp(-2*k3*freq))

## 7. Return State Space Form

ssf.mod = list(mDelta=mDelta, mPhi=mPhi, mOmega=mOmega)
CheckSsf(ssf.mod)

}

start = c(0.1,0.1,0.1,0.3,0.3,0.3,0.06,0.005,0.005,0.005,0.001,0.001,0.001,
0.003,0.1,0.003,0.01)
names(start)=c("logkappa1","logkappa2","logkappa3","lambda1","lambda2","lambda3","logDelta","logSigma1","logSigma2","logSigma3","Corr12","Corr13","Corr23","logSigmat1","logSigmat2","logSigmat3","logSigmat4")

tau = c(0.25,0.5,1,5)

l = c(0,0,0,0,0,0,0,0,0,-0,-1,-1,-1,0,0,0,0)
u = c(10,10,10,10,10,10,0.5,10,10,10,1,1,1,10,10,10,10)


ans.vasicek = SsfFit(start, fama.bliss, genvasicek.ssf, tau=tau, freq=1/12,
trace=T, lower=l, upper=u, control=nlminb.control(abs.tol=1e-6,,
rel.tol=1e-6, x.tol=1e-6, eval.max=1000,iter.max=135))


On 8/13/07, Bernardo Ribeiro <[EMAIL PROTECTED]> wrote:
>
> Hey all,
>
> I am trying to work under a State Space form, but I didn't get the help
> exactly.
> Have anyone eles used this functions?
>
> I was used to work with S-PLUS, but I have some codes I need to adpt.
>
> Thanks alot,
>
> Bernardo
>

[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Using sunflowerplot to add points in a xyplot panel

[Ronaldo Reis Junior]

Hi,

I use panel.points to add points to a xyplot graphic. But I like to use the 
sunflowerplot to plot my points because this is very superimposed. It is 
possible to use this? I try but it dont work directly. It may be need to put 
this function inside a panel.??? 

Thanks
Ronaldo
-- 
Where there's a will, there's a relative.
--
> Prof. Ronaldo Reis Júnior
|  .''`. UNIMONTES/Depto. Biologia Geral/Lab. de Ecologia
| : :'  : Campus Universitário Prof. Darcy Ribeiro, Vila Mauricéia
| `. `'` CP: 126, CEP: 39401-089, Montes Claros - MG - Brasil
|   `- Fone: (38) 3229-8187 | [EMAIL PROTECTED] | [EMAIL PROTECTED]
| http://www.ppgcb.unimontes.br/ | ICQ#: 5692561 | LinuxUser#: 205366

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] invert 160000x160000 matrix

[Jiao Yang]

Tir,

Thank you very much for the note.  I'm using an algorithm to analyze a data set 
of 400 variables and the algorithm uses the inverse of a 16x16 positive 
definite matrix.  The matrix is stored as a text file.

Can you please give some reference to  in-place inverse technique?  I googled 
it, yet did not get any good results.  Thank you for your time.

Best regards,
Jiao


- Original Message -
From: "Patnaik, Tirthankar " <[EMAIL PROTECTED]>
Date: Tuesday, August 14, 2007 2:07 am
Subject: RE: [R] invert 16x16 matrix
To: Moshe Olshansky <[EMAIL PROTECTED]>, Paul Gilbert <[EMAIL PROTECTED]>, Jiao 
Yang <[EMAIL PROTECTED]>
Cc: r-help@stat.math.ethz.ch


> A variety of tricks would need to be used to invert a matrix of this 
> size. If there are any other properties of the matrix that you know 
> (symmetric, positive definite, etc, sparse) then they could be useful 
> too. You could partition the matrix first, then use an in-place 
> inverse technique for each block to individually calculate the 
> blocks-inverses, then combine to get the inverse of the initial 
> matrix. Again, if the implementation is actually solving an Ax-B = 0 
> system of equations, then there are specific methods for these too, 
> like an LU decomp, for instance. You might also want to check out some 
> texts for this, like the Numerical Recipes. 
> How's the matrix stored right now?
> 
> Best,
> -Tir
> 
> 
> 
> Tirthankar Patnaik
> India Strategy
> Citigroup Investment Research
> +91-22-6631 9887
> 
> For important disclosures regarding Citigroup Investment Research, 
> including with respect to any issuers mentioned herein, please refer 
> to the Citigroup Investment Research disclosure website at
> 
> 
> > -Original Message-
> > From: [EMAIL PROTECTED] 
> > [ On Behalf Of Moshe Olshansky
> > Sent: Tuesday, August 14, 2007 6:40 AM
> > To: Paul Gilbert; Jiao Yang
> > Cc: r-help@stat.math.ethz.ch
> > Subject: Re: [R] invert 16x16 matrix
> > 
> > While inverting the matrix may be a problem, if you need to 
> > solve an equation A*x = b you do not need to invert A, there 
> > exist iterative methods which do need A or inv(A) - all you 
> > need to provide is a function that computes A*x for an 
> > arbitrary vector x.
> > For such a large matrix this may be slow but possible.
> > 
> > --- Paul Gilbert <[EMAIL PROTECTED]>
> > wrote:
> > 
> > > I don't think you can define a matrix this large in R, even if you 
> 
> > > have the memory. Then, of course, inverting it there may be other 
> 
> > > programs that have limitations.
> > > 
> > > Paul
> > > 
> > > Jiao Yang wrote:
> > > > Can R invert a 16x16 matrix with all
> > > positive numbers?  Thanks a lot!
> > > > 
> > > > __
> > > > R-help@stat.math.ethz.ch mailing list 
> > > > 
> > > > PLEASE do read the posting guide
> > > 
> > > > and provide commented, minimal, self-contained,
> > > reproducible code.
> > >
> > ==
> > ==
> > > 
> > > La version française suit le texte anglais.
> > > 
> > >
> > --
> > --
> > > 
> > > This email may contain privileged and/or confidential 
> > > inform...{{dropped}}
> > > 
> > > __
> > > R-help@stat.math.ethz.ch mailing list
> > > 
> > > PLEASE do read the posting guide
> > > 
> > > and provide commented, minimal, self-contained, reproducible code.
> > >
> > 
> > __
> > R-help@stat.math.ethz.ch mailing list
> > 
> > PLEASE do read the posting guide 
> > 
> > and provide commented, minimal, self-contained, reproducible code.
> >

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] labelling plots with ancillary data in data.frame

[Wesley Roberts]

Hi All,

I am busy using R to do some regression modelling and have been using 
plot(x,y,"") to visualise my variables. I would now like to label my points 
using data stored in the data.frame used for the regression analysis. For 
example each of my data points is made up of a field measured forest volume 
value and a remotely sensed vegetation estimate (NDVI). Each point is an 
enumeration plot and I would like to label each the points in the 
xy-scatterplot with their respective plot numbers. Is this possible in R, if so 
how do I go about doing it?

Many thanks for your help

Wesley

Wesley Roberts MSc.
Researcher: Forest Assessment (Remote Sensing & GIS)
Forestry and Forest Products Research Centre
CSIR
Tel: +27 (31) 242-2353
Fax: +27 (31) 261-1216
http://ffp.csir.co.za/

"To know the road ahead, ask those coming back."
- Chinese proverb


-- 
This message is subject to the CSIR's copyright, terms and conditions and
e-mail legal notice. Views expressed herein do not necessarily represent the
views of the CSIR.
 
CSIR E-mail Legal Notice
http://mail.csir.co.za/CSIR_eMail_Legal_Notice.html 
 
CSIR Copyright, Terms and Conditions
http://mail.csir.co.za/CSIR_Copyright.html 
 
For electronic copies of the CSIR Copyright, Terms and Conditions and the CSIR
Legal Notice send a blank message with REQUEST LEGAL in the ...{{dropped}}

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Problem with "by": does not work with ttest (but with lme)

[Chuck Cleland]

Daniel Stahl wrote:
> Hello,
> 
> I would like to do a large number of e.g. 1000 paired ttest using the 
> by-function. But instead of using only the data within the 1000 groups, R 
> caclulates 1000 times the ttest for the full data set(The same happens with 
> Wilcoxon test). However, the by-function works fine with the lme function.
> Did I just miss something or is it really not working? If not, is there any 
> other possibility to avoid loops? 
> Thanks 
> Daniel
> 
> Here is the R help example for "by" 
>  require(stats)
>  attach(warpbreaks)
>  by(warpbreaks, tension, function(x) lm(breaks ~ wool, data = x))
> *->works great
> by(warpbreaks,tension,function(x)t.test(breaks ~ wool,data=warpbreaks,paired 
> = TRUE))
> *Same output for each level of tension:
> 
> tension: L
> 
>   Paired t-test
> 
> data:  breaks by wool
> t = 1.9956, df = 26, p-value = 0.05656
> alternative hypothesis: true difference in means is not equal to 0
> 95 percent confidence interval:
> -0.1735803 11.7291358
> sample estimates:
> mean of the differences
> 5.78
> 
> 
> 
> tension: M
> 
>   Paired t-test
> 
> data:  breaks by wool
> t = 1.9956, df = 26, p-value = 0.05656
> alternative hypothesis: true difference in means is not equal to 0
> 95 percent confidence interval:
> -0.1735803 11.7291358
> sample estimates:
> mean of the differences
> 5.78
> 
> 
> 
> tension: H
> 
>   Paired t-test
> 
> data:  breaks by wool
> t = 1.9956, df = 26, p-value = 0.05656
> alternative hypothesis: true difference in means is not equal to 0
> 95 percent confidence interval:
> -0.1735803 11.7291358
> sample estimates:
> mean of the differences
> 5.78

  Try something like this:

library(MASS)
df <- mvrnorm(30, mu=c(-1,1), Sigma = diag(2))
df <- as.data.frame(df)
df$GROUP <- rep(1:3, each=10)

df.uni <- reshape(df, varying = list(c("V1","V2")), v.names="Y",
direction="long")

by(df.uni, df.uni$GROUP, function(x)t.test(Y ~ time,
   data = x, paired = TRUE))

df.uni$GROUP: 1

Paired t-test

data:  Y by time
t = -4.3719, df = 9, p-value = 0.001792
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -3.249894 -1.033530
sample estimates:
mean of the differences
  -2.141712

---

df.uni$GROUP: 2

Paired t-test

data:  Y by time
t = -6.4125, df = 9, p-value = 0.0001234
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -3.277425 -1.568074
sample estimates:
mean of the differences
  -2.422749

---

df.uni$GROUP: 3

Paired t-test

data:  Y by time
t = -4.4918, df = 9, p-value = 0.001507
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -3.581428 -1.182313
sample estimates:
mean of the differences
  -2.381871

> --
> 
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code. 

-- 
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Question about unicode characters in tcltk

[R Help]

hello list,

Can someone help me figure out why the following code doesn't work?
I'm trying to but both Greek letters and subscripts into a tcltk menu.
 The code creates all the mu's, and the 1 and 2 subscripts, but it
won't create the 0.  Is there a certain set of characters that R won't
recognize the unicode for?  Or am I input the \u2080 incorrectly?

library(tcltk)
m <-tktoplevel()
frame1 <- tkframe(m)
frame2 <- tkframe(m)
frame3 <- tkframe(m)
entry1 <- tkentry(frame1,width=5,bg='white')
entry2 <- tkentry(frame2,width=5,bg='white')
entry3 <- tkentry(frame3,width=5,bg='white')

tkpack(tklabel(frame1,text='\u03bc\u2080'),side='left')
tkpack(tklabel(frame2,text='\u03bc\u2081'),side='left')
tkpack(tklabel(frame3,text='\u03bc\u2082'),side='left')

tkpack(frame1,entry1,side='top')
tkpack(frame2,entry2,side='top')
tkpack(frame3,entry3,side='top')

thanks
-- Sam

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Error message when using zero-inflated count regression model in package zicounts

[Hank Stevens]

Hi Jim,
Two thoughts come to me, unencumbered by the thought process or  
knowledge of zicounts:
1. Is Poisson really NOT appropriate? (do you have to use zicounts?)
2. Are you 110% certain that all variables are the same length? Would  
NA's interfere?
Cheers,
Hank
On Aug 13, 2007, at 5:10 PM, James Milks wrote:

> I have data on number of vines per tree for ~550 trees.  Over half of
> the trees did not have any vines and the data is fairly skewed
> (median = 0, mean = 1.158, 3rd qu. = 1.000).  I am attempting to
> investigate whether plot location (four sites), species (I'm using
> only the four most common species), or tree dbh has a significant
> influence on the number of vines per tree.  When I attempted to use
> the zicounts function, R gave me the following error message:
>
>> vines.zip<-zicounts(resp=Total.vines~.,x=~Site+Species+DBH,z=~Site
> +Species+DBH,distrname="ZIP",data=sycamores.1)
> Error in ifelse(y == 0, 1, y/mu) : dim<- : dims [product 12] do not
> match the length of object [549]
> In addition: Warning messages:
> 1: longer object length
> is not a multiple of shorter object length in: x[good, ] * w
> 2: longer object length
> is not a multiple of shorter object length in: eta + offset
> 3: longer object length
> is not a multiple of shorter object length in: y/mu
>
> I do not know enough about the calculations done in the function to
> interpret the error messages.  Is there a glitch in my data and if
> yes, what is it?
>
> Thanks for your help.
>
> Jim Milks
>
> Graduate Student
> Environmental Sciences Ph.D. Program
> 136 Biological Sciences
> Wright State University
> 3640 Colonel Glenn Hwy
> Dayton, OH 45435
>
>
>
> [[alternative HTML version deleted]]
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting- 
> guide.html
> and provide commented, minimal, self-contained, reproducible code.



Dr. Hank Stevens, Associate Professor
338 Pearson Hall
Botany Department
Miami University
Oxford, OH 45056

Office: (513) 529-4206
Lab: (513) 529-4262
FAX: (513) 529-4243
http://www.cas.muohio.edu/~stevenmh/
http://www.muohio.edu/ecology/
http://www.muohio.edu/botany/

"E Pluribus Unum"

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Odp: Problem with "by": does not work with ttest (but with lme)

[Petr PIKAL]

Hi

[EMAIL PROTECTED] napsal dne 14.08.2007 15:11:11:

> Hello,
> 
> I would like to do a large number of e.g. 1000 paired ttest using the 
by-
> function. But instead of using only the data within the 1000 groups, R 
> caclulates 1000 times the ttest for the full data set(The same happens 
with 
> Wilcoxon test). However, the by-function works fine with the lme 
function.
> Did I just miss something or is it really not working? If not, is there 
any 
> other possibility to avoid loops? 
> Thanks 
> Daniel
> 
> Here is the R help example for "by" 
>  require(stats)
>  attach(warpbreaks)
>  by(warpbreaks, tension, function(x) lm(breaks ~ wool, data = x))
> *->works great
> by(warpbreaks,tension,function(x)t.test(breaks ~ 
wool,data=warpbreaks,paired = TRUE))

What about

by(warpbreaks,tension,function(x)t.test(breaks ~ wool,data=x,paired = 
TRUE))

Regards
Petr


> *Same output for each level of tension:
> 
> tension: L
> 
>Paired t-test
> 
> data:  breaks by wool
> t = 1.9956, df = 26, p-value = 0.05656
> alternative hypothesis: true difference in means is not equal to 0
> 95 percent confidence interval:
> -0.1735803 11.7291358
> sample estimates:
> mean of the differences
> 5.78
> 
> 
> 
> tension: M
> 
>Paired t-test
> 
> data:  breaks by wool
> t = 1.9956, df = 26, p-value = 0.05656
> alternative hypothesis: true difference in means is not equal to 0
> 95 percent confidence interval:
> -0.1735803 11.7291358
> sample estimates:
> mean of the differences
> 5.78
> 
> 
> 
> tension: H
> 
>Paired t-test
> 
> data:  breaks by wool
> t = 1.9956, df = 26, p-value = 0.05656
> alternative hypothesis: true difference in means is not equal to 0
> 95 percent confidence interval:
> -0.1735803 11.7291358
> sample estimates:
> mean of the differences
> 5.78
> 
> 
> 
> 
> 
> 
> 
> --
> 
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide 
http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Problem with "by": does not work with ttest (but with lme)

[Daniel Stahl]

Hello,

I would like to do a large number of e.g. 1000 paired ttest using the 
by-function. But instead of using only the data within the 1000 groups, R 
caclulates 1000 times the ttest for the full data set(The same happens with 
Wilcoxon test). However, the by-function works fine with the lme function.
Did I just miss something or is it really not working? If not, is there any 
other possibility to avoid loops? 
Thanks 
Daniel

Here is the R help example for "by" 
 require(stats)
 attach(warpbreaks)
 by(warpbreaks, tension, function(x) lm(breaks ~ wool, data = x))
*->works great
by(warpbreaks,tension,function(x)t.test(breaks ~ wool,data=warpbreaks,paired = 
TRUE))
*Same output for each level of tension:

tension: L

Paired t-test

data:  breaks by wool
t = 1.9956, df = 26, p-value = 0.05656
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
-0.1735803 11.7291358
sample estimates:
mean of the differences
5.78



tension: M

Paired t-test

data:  breaks by wool
t = 1.9956, df = 26, p-value = 0.05656
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
-0.1735803 11.7291358
sample estimates:
mean of the differences
5.78



tension: H

Paired t-test

data:  breaks by wool
t = 1.9956, df = 26, p-value = 0.05656
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
-0.1735803 11.7291358
sample estimates:
mean of the differences
5.78







--

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] papply, reuse of data

[Markus Schmidberger]

Hello,

we are doing some computations with MPI and Rmpi. Therefore papply is a 
very comfortable function.
Is it possible to reuse calculated data on slaves?

We do some calculation with papply,
then we have to do some overall calculation at the master
then we want to do some calculation on the slaves with the old results 
and some new results from the master.

We do not want to send all data again from master to slave (time 
expensive). Is there a way to do this with papply?

Thanks
Markus

-- 
Dipl.-Tech. Math. Markus Schmidberger

Ludwig-Maximilians-Universität München
IBE - Institut für medizinische Informationsverarbeitung,
Biometrie und Epidemiologie
Marchioninistr. 15, D-81377 Muenchen
URL: http://ibe.web.med.uni-muenchen.de 
Mail: Markus.Schmidberger [at] ibe.med.uni-muenchen.de

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] graph dimensions default

[Prof Brian Ripley]

On Tue, 14 Aug 2007, Simon Pickett wrote:

> Yes,
>
> Thankyou, that does the trick nicely. I thought that kind of thing could
> be specified using par() but I guess not.

As I said, size is not a property of the plot.
And par() applies to the current device, not future ones.

>
> Thanks again.
>
>
>
>> On Tue, 14 Aug 2007, Simon Pickett wrote:
>>
>>> Hi,
>>>
>>> I would like to (if possible) set the default width and height for
>>> graphs
>>> at the start of each session and have each new graphic device overwrite
>>> the previous one.
>>
>> Hmm.  It is graphics devices that have dimensions, and plots that
>> overwrite other plots on a device, so your intentions are pretty unclear.
>> (If you resize a device window the plot dimensions change so they are not
>> intrinsic to the plot.)
>>
>> If you want the default behaviour to be like normal but with, say, a wider
>> onscreen device window you can have (on Windows, which you didn't say)
>>
>> mywindows <- function(...) windows(width=10, height=6, ...)
>> options(device="mywindows")
>>
>> in your ~/.Rprofile .  Otherwise, please try again to tell us what you
>> do want.
>>
>>
>>>
>>> I only know how to do this using windows(width=,height=...) which opens
>>> up
>>> a new plotting device every time, so I end up with lots of graphs all
>>> over
>>> the place until I get the one I want!
>>>
>>> Thanks in advance,
>>>
>>> Simon
>>>
>>>
>>> Simon Pickett
>>> PhD student
>>> Centre For Ecology and Conservation
>>> Tremough Campus
>>> University of Exeter in Cornwall
>>> TR109EZ
>>> Tel 01326371852
>>>
>>> http://www.uec.ac.uk/biology/research/phd-students/simon_pickett.shtml
>>>
>>> __
>>> R-help@stat.math.ethz.ch mailing list
>>> https://stat.ethz.ch/mailman/listinfo/r-help
>>> PLEASE do read the posting guide
>>> http://www.R-project.org/posting-guide.html
>>> and provide commented, minimal, self-contained, reproducible code.
>>>
>>
>> --
>> Brian D. Ripley,  [EMAIL PROTECTED]
>> Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
>> University of Oxford, Tel:  +44 1865 272861 (self)
>> 1 South Parks Road, +44 1865 272866 (PA)
>> Oxford OX1 3TG, UKFax:  +44 1865 272595
>>
>
>
> Simon Pickett
> PhD student
> Centre For Ecology and Conservation
> Tremough Campus
> University of Exeter in Cornwall
> TR109EZ
> Tel 01326371852
>
> http://www.uec.ac.uk/biology/research/phd-students/simon_pickett.shtml
>
>

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Linear Regression with slope equals 0

[Ted Harding]

On 14-Aug-07 11:36:37, [EMAIL PROTECTED] wrote:
> 
>  Hi there, am trying to run a linear regression with a slope of 0.
> 
>  I have a dataset as follows
> 
>  t d
>  1 303
>  2 302
>  3 304
>  4 306
>  5 307
>  6 303
> 
>  I would like to test the significance that these points would lie on a
> horizontal straight line.
> 
>  The standard regression lm(d~t) doesn't seem to allow the slope to be
> set.

The model d~1 will fit a constant (the mean), i.e. a regressio with
slope = 0. The model d~t will fit the usual linear regression.
The two con be compared with anova(), as well as getting the details
of the individual fits with summary().

E.g. (with your example):

  d<-c(303,302,304,306,207,303)
  t<-c(1,2,3,4,5,6)
  lm0<-lm(u~1);lm1<-lm(u~t);anova(lm0,lm1)

##Analysis of Variance Table
##Model 1: u ~ 1
##Model 2: u ~ t
##  Res.DfRSS Df Sum of Sq  F Pr(>F)
##1  5 7785.5   
##2  4 6641.4  11144.1 0.6891 0.4531

summary(lm0)
## Call: lm(formula = u ~ 1)
## Residuals:
##1 2 3 4 5 6 
## 15.5  14.5  16.5  18.5 -80.5  15.5 
## Coefficients:
##Estimate Std. Error t value Pr(>|t|)
## (Intercept)   287.50  16.11   17.85 1.01e-05 ***
##Residual standard error: 39.46 on 5 degrees of freedom

mean(d)
## [1] 287.5

summary(lm1)
## Call: lm(formula = u ~ t)
## Residuals:
##  1   2   3   4   5   6 
## -4.714   2.371  12.457  22.543 -68.371  35.714 
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)   
## (Intercept)  315.800 37.934   8.325  0.00114 **
## t -8.086  9.740  -0.830  0.45314   
## Residual standard error: 40.75 on 4 degrees of freedom
## Multiple R-Squared: 0.147,  Adjusted R-squared: -0.0663 
## F-statistic: 0.6891 on 1 and 4 DF,  p-value: 0.4531 

The P-value for the slope in lm1 is the same as the P-value
returned by anova().

If you want to force a particular non-zero slope (e.g. s0)
for comparison with the data, you can use

  lm0 <- lm(d - s0*t ~ 1),

compared with

  lm1<- lm(d- s0*t ~ t)

for instance.

Hoping this helps,
Ted.


E-Mail: (Ted Harding) <[EMAIL PROTECTED]>
Fax-to-email: +44 (0)870 094 0861
Date: 14-Aug-07   Time: 13:16:05
-- XFMail --

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] vertically oriented color key in heatmaps

[Rajarshi Guha]

On Aug 14, 2007, at 5:11 AM, Jim Lemon wrote:

> Rajarshi Guha wrote:
>> Hi, I have some data which I was plotting using image(). I wanted  
>> to  add a vertical color key to the plot and I found that heatmap. 
>> 2 in  gplots does let me add a color key. However, I was thinking  
>> of a  vertical bar with the color range rather than  the style  
>> that gplots  provides.
>> Is there any package (or code snippet) that would let me add a   
>> vertical color key to an image() or heatmap plot?
> Hi Rajarshi,
> Have a look at color.legend in the plotrix package.

Thanks!

---
Rajarshi Guha  <[EMAIL PROTECTED]>
GPG Fingerprint: 0CCA 8EE2 2EEB 25E2 AB04  06F7 1BB9 E634 9B87 56EE
---
Artificial intelligence has the same relation to intelligence as
artificial flowers have to flowers.
-- David Parnas

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] graph dimensions default

[Simon Pickett]

Yes,

Thankyou, that does the trick nicely. I thought that kind of thing could
be specified using par() but I guess not.

Thanks again.



> On Tue, 14 Aug 2007, Simon Pickett wrote:
>
>> Hi,
>>
>> I would like to (if possible) set the default width and height for
>> graphs
>> at the start of each session and have each new graphic device overwrite
>> the previous one.
>
> Hmm.  It is graphics devices that have dimensions, and plots that
> overwrite other plots on a device, so your intentions are pretty unclear.
> (If you resize a device window the plot dimensions change so they are not
> intrinsic to the plot.)
>
> If you want the default behaviour to be like normal but with, say, a wider
> onscreen device window you can have (on Windows, which you didn't say)
>
> mywindows <- function(...) windows(width=10, height=6, ...)
> options(device="mywindows")
>
> in your ~/.Rprofile .  Otherwise, please try again to tell us what you
> do want.
>
>
>>
>> I only know how to do this using windows(width=,height=...) which opens
>> up
>> a new plotting device every time, so I end up with lots of graphs all
>> over
>> the place until I get the one I want!
>>
>> Thanks in advance,
>>
>> Simon
>>
>>
>> Simon Pickett
>> PhD student
>> Centre For Ecology and Conservation
>> Tremough Campus
>> University of Exeter in Cornwall
>> TR109EZ
>> Tel 01326371852
>>
>> http://www.uec.ac.uk/biology/research/phd-students/simon_pickett.shtml
>>
>> __
>> R-help@stat.math.ethz.ch mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide
>> http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>>
>
> --
> Brian D. Ripley,  [EMAIL PROTECTED]
> Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
> University of Oxford, Tel:  +44 1865 272861 (self)
> 1 South Parks Road, +44 1865 272866 (PA)
> Oxford OX1 3TG, UKFax:  +44 1865 272595
>


Simon Pickett
PhD student
Centre For Ecology and Conservation
Tremough Campus
University of Exeter in Cornwall
TR109EZ
Tel 01326371852

http://www.uec.ac.uk/biology/research/phd-students/simon_pickett.shtml

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Odp: Linear Regression with slope equals 0

[Petr PIKAL]

Hi

Not being a trained statistician regression with slope = 0 seems odd to 
me.

If you do

> fit<-lm(d~t)
> summary(fit)

Call:
lm(formula = d ~ t)

Residuals:
   123456 
 0.04762 -1.43810  0.07619  1.59048  2.10476 -2.38095 

Coefficients:
Estimate Std. Error t value Pr(>|t|) 
(Intercept) 302.4667 1.7849  169.45 7.28e-09 ***
t 0.4857 0.45831.060.349 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Residual standard error: 1.917 on 4 degrees of freedom
Multiple R-Squared: 0.2192, Adjusted R-squared: 0.02402 
F-statistic: 1.123 on 1 and 4 DF,  p-value: 0.349 

you will get estimate of your coeficients and AFAIK they are tested 
against Ho that they differ from 0.

So as you see your "t" is not statistically different from 0.

Regards
Petr


[EMAIL PROTECTED] napsal dne 14.08.2007 13:36:37:

> 
>  Hi there, am trying to run a linear regression with a slope of 0.
> 
>  I have a dataset as follows
> 
>  t d
>  1 303
>  2 302
>  3 304
>  4 306
>  5 307
>  6 303
> 
>  I would like to test the significance that these points would lie on a
> horizontal straight line.
> 
>  The standard regression lm(d~t) doesn't seem to allow the slope to be 
set.
> 
>  Any help very welcome.
> 
>  ed
> 
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide 
http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Import of Access data via RODBC changes column name ("NO" to "Expr1014") and the content of the column

[Maciej Hoffman-Wecker]

Dear Professor Ripley,

Thank you very much for your response. I send the problem, as I didn't have any 
more ideas were to search for the reason. I didn't say this is a R bug, knowing 
the responses on such mails.-) 

But I succeeded in developing a tiny example, that reproduces the bug (wherever 
it is).

I generated a small Access data base "test2.mdb" with one table "Tab1" and 
following columns:

"Field name""Field type"
  F1Number
  NONumber
  F2Number

(sorry if the Access identifiers are not the standard ones, as I have a german 
Access version)

The content of the "Tab1" table is:

F1  NO  F2
1   1   1
2   2   2
0   1
1   0   0

(The column "NO" contains one missing)

Now if I import the data into R, I get the following results:

> library(RODBC)
> .con <- odbcConnectAccess("./test2.mdb")
> (.d <- try(sqlQuery(.con, "select * from Tab1")))
  F1 NO F2
1  1  1  1
2  2  2  2
3  0 NA  1
4  1  0  0
> (.d <- try(sqlQuery(.con, "select F1 , NO , F2 from Tab1")))
  F1 Expr1001 F2
1  10  1
2  20  2
3  00  1
4  10  0
> close(.con)

So the problem occurs if the column names are specified within the query.
Is the query "select F1 , NO , F2 from Tab1" invalid?

Regarding the memory issue, I _knew_ that there must be a reason for the 
running out of memory space. Sorry for not being more specific. My question 
than is:

Is there a way to 'reset' the environment without quitting R and restarting it?

Thank you for your help. 

Kind regards,
Maciej


-Ursprüngliche Nachricht-
Von: Prof Brian Ripley [mailto:[EMAIL PROTECTED] 
Gesendet: Dienstag, 14. August 2007 11:51
An: Maciej Hoffman-Wecker
Cc: r-help@stat.math.ethz.ch
Betreff: Re: [R] Import of Access data via RODBC changes column name ("NO" to 
"Expr1014") and the content of the column

On Tue, 14 Aug 2007, Maciej Hoffman-Wecker wrote:

>
> Dear all,
>
> I have some problems with importing data from an Access data base via 
> RODBC to R. The data base contains several tables, which all are 
> imported consecutively. One table has a column with column name "NO". 
> If I run the code attached on the bottom of the mail I get no 
> complain, but the column name (name of the respective vector of the 
> data.frame) is "Expr1014" instead of "NO". Additionally the original 
> column (type
> "text") containes "0"s and missings, but the imported column contains 
> "0"s only (type "int"). If I change the column name in the Access data 
> base to "NOx", the import works fine with the right name and the same 
> data.
>
> Previously I generated a tiny Access data base which reproduced the 
> problem. To be on the safe site I installed the latest version (2.5.1) 
> and now the example works fine, but within my production process the 
> error still remaines. An import into excel via ODBC works fine.
>
> So there is no way to figure it out whether this is a bug or a
> feature.-)

It's most likely an ODBC issue, but you have not provided a reproducible 
example.

> The second problem I have is that when I rerun "rm(list = ls(all = 
> T)); gc()" and the import several times I get the following error:
>
> Error in odbcTables(channel) : Calloc could not allocate (263168 of 1) 
> memory In addition: Warning messages:
> 1: Reached total allocation of 447Mb: see help(memory.size) in:
> odbcQuery(channel, query, rows_at_time)
> 2: Reached total allocation of 447Mb: see help(memory.size) in:
> odbcQuery(channel, query, rows_at_time)
> 3: Reached total allocation of 447Mb: see help(memory.size) in:
> odbcTables(channel)
> 4: Reached total allocation of 447Mb: see help(memory.size) in:
> odbcTables(channel)
>
> which is surprising to me, as the first two statements should delete 
> all

How do you _know _what they 'should' do?  That only deletes all objects in the 
workspace, not all objects in R, and not all memory blocks used by R.

Please do read ?"Memory-limits" for the possible reasons.

Where did '447Mb' come from?  If this machine has less than 2Gb of RAM, buy 
some more.


> objects and recover the memory. Is this only a matter of memory? Is
> there any logging that reduces the memory? Or is this issue connected to
> the upper problem?
>
> I added the code on the bottom - maybe there is some kind of misuse I
> lost sight of. Any hints are appreciated.
>
> Kind regards,
> Maciej
>
>> version
>   _
> platform   i386-pc-mingw32
> arch   i386
> os mingw32
> system i386, mingw32
> status
> major  2
> minor  5.1
> year   2007
> month  06
> day27
> svn rev42083
> language   R
> version.string R version 2.5.1 (2007-06-27)
>
>
> ## code
>
> get.table <- function(name, db, drop = NULL){
>  .con <- try(odbcConnectAccess(db), silent = T)
>  if(!inherits(.con, "RODBC")) return(.con)
>  ## exclude memo columns
>  .t <- try(sqlColumns(.con, name))
>  

Re: [R] Linear Regression with slope equals 0

[Prof Brian Ripley]

On Tue, 14 Aug 2007, [EMAIL PROTECTED] wrote:

>
> Hi there, am trying to run a linear regression with a slope of 0.
>
> I have a dataset as follows
>
> t d
> 1 303
> 2 302
> 3 304
> 4 306
> 5 307
> 6 303
>
> I would like to test the significance that these points would lie on a
> horizontal straight line.
>
> The standard regression lm(d~t) doesn't seem to allow the slope to be set.

lm(d ~ 1) does, though, to zero.

More generally you can use offset(), e.g. lm(d ~ offset(7*t)) forces a 
slope of 7.

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] graph dimensions default

[Prof Brian Ripley]

On Tue, 14 Aug 2007, Simon Pickett wrote:

> Hi,
>
> I would like to (if possible) set the default width and height for graphs
> at the start of each session and have each new graphic device overwrite
> the previous one.

Hmm.  It is graphics devices that have dimensions, and plots that 
overwrite other plots on a device, so your intentions are pretty unclear. 
(If you resize a device window the plot dimensions change so they are not 
intrinsic to the plot.)

If you want the default behaviour to be like normal but with, say, a wider 
onscreen device window you can have (on Windows, which you didn't say)

mywindows <- function(...) windows(width=10, height=6, ...)
options(device="mywindows")

in your ~/.Rprofile .  Otherwise, please try again to tell us what you 
do want.


>
> I only know how to do this using windows(width=,height=...) which opens up
> a new plotting device every time, so I end up with lots of graphs all over
> the place until I get the one I want!
>
> Thanks in advance,
>
> Simon
>
>
> Simon Pickett
> PhD student
> Centre For Ecology and Conservation
> Tremough Campus
> University of Exeter in Cornwall
> TR109EZ
> Tel 01326371852
>
> http://www.uec.ac.uk/biology/research/phd-students/simon_pickett.shtml
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] makeSOCKcluster

[Luke Tierney]

Try using the option homogenoeus=FALSE and make sure the appropriate
environment variables are set on the worker nodes.

Best,

luke

On Tue, 14 Aug 2007, Michael Janis wrote:

> Hi,
>
>
>
> I am attempting to implement a mixed (windows/linux) snow sockets
> parallelism in R, but am running into difficulties similar to a post made
> Aug 31, 2006 under the same subject heading.  I feel like I may be one or
> two non-obvious steps away from getting it all working, but I'm stuck.  If
> anyone can shed some light on this (I believe Prof. Tierney stated that he
> has successfully run a snow master on Windows with slave nodes on Linux
> using ssh.exe through Cygwin, which is exactly what I am attempting), I'd be
> most grateful.
>
>
>
> SNOW master on Windows information:
>
> WinXP 32bit
>
> R-2.5.1
>
> I built a windows package (.zip) version of snow-0.2-3 for R using Rtools
> and installed it without trouble.
>
> Cygwin including ssh as an exe and in the path
>
>
>
> Linux slave nodes information:
>
> ROCKS compute nodes, each 64bit (I mention this difference between the win
> and linux platforms out of desperation but I don't think it's an issue since
> snow should be agnostic).  (btw, the SOCK version of snow here is just for
> me to get used to parallelism in R before attempting MPI).
>
> R-2.5.1, local directory install (no su, global R install is older).
>
> Snow installed to R-2.5.1 instance
>
> Unattended ssh authentication through public/private key pairs.  Each node
> allows since they are all NFS.
>
>
>
> Rgui commands:
>
>> library(snow)
>
>> cl<-makeCluster(c("ip.path.to.node1","ip.path.to.node2"),type="SOCK")
>
> Result:
>
> R becomes unresponsive and has to be forcibly closed.  Sometimes before that
> happens message stating env: C:\pathToR\library\snow\RSOCKnode.sh no file
> found appears (if I make R redraw the interface screen)
>
> Troubleshooting step1:
>
> Interspersed print commands within C:\localRinstall\library\snow\R\snow
> file.
>
> Result of troubleshooting step1:
>
>It appears that this snow file hangs on line:
> system(paste(rshcmd, "-l", user, machine, "env", env, script))
>
> Resolution of troubleshooting step1:
>
>  Manually attempted paste commands in R:
>
>  >system("ssh -l me ip.path.to.node1 ls")
>
>Result: gets ls of remote machine; so the script path must be
> incorrect and the offending command
>
>  Located RSOCKnode.sh on remote machine (for all remote machines, since
> it's NFS)
>
>Inserted line into C:\localRinstall\library\snow\R\snow file
> just before the offending line system(paste(rshcmd, "-l", user, machine,
> "env", env, script)).  That inserted line hard-codes the path of the script
> and the script name like so:
> script<-"/rootpath/home/cluster/me/R-2.5.1/library/snow/RSOCKnode.sh".  This
> worked, and the next print statement following the system call in the snow
> file now prints to the screen.  But again R hangs, this time at the line:
> con <- socketConnection(port = port, server=TRUE, blocking=TRUE,
> open="a+b").
>
> Troubleshooting step2:
>
>  Attempted manual socketConnection via R commands:
>
>  > con<-socketConnection("ip.path.to.node1",port=22)
>
>> showConnections()
>
>description classmode text   isopen   can read can
> write
>
> 3  "->ip.path.to.node1" "socket" "a+" "text" "opened" "yes""yes"
>
> *** so the socketConnection does work, but in SNOW it hangs at this last
> command, and I'm completely stuck.
>
>
>
> As a parallel test, I have installed SNOW on the R-2.5.1 instance on the
> linux cluster.  If I access the head node and launch R, I can use SNOW with
> sockets successfully (from the linux head node as master to the linux
> compute nodes as slaves). However, running R more or less continually on the
> head node of a cluster is bad form.  Ideally I would like to run a windows
> snow master to the linux node slaves, but I can't seem to get past this
> point.
>
>
>
> For reference, I include the original post below.  I am stuck at one command
> in the snow file past where this last post had problems.
>
>
>
> I hope this is clear.  I am intent on solving this problem, so feel free to
> ask questions if you have feedback but my description is not clear.  I
> really appreciate any help!
>
>
>
> Yours,
>
>
>
> Michael Janis
>
> UCLA Bioinformatics
>
>
>
> *transcript from original post entitles "makeSOCKcluster" follows*
>
>
>
>
>
> makeSOCKcluster
>
> Click to flag this post
>
>
>
> by Hrishikesh Rajpathak Aug 31, 2006; 10:39pm :: Rate this Message: - Use
> ratings to moderate (?)
>
>
>
> Reply | Reply to Author | View Threaded | Show Only this Message
>
> Hi,
>
>
>
>  I am a newbie to R and trying to implement parallelism in R. I am
> currently using R-2.3.1, and Cygwin to run R on Windows xp.
>
>
>
>  ssh and all are working fine,
>
>
>
>  When I try to create a socket connection as
>
>
>
>   makeSOCKcluster(c("localhost","localhost")),
>
>
>
>  it ju

Re: [R] diffusing GIS data in maps

[Felix Andrews]

Package akima allows simple 2D interpolation from points to a grid.


On 8/14/07, Lawrence D. Brenninkmeyer <[EMAIL PROTECTED]> wrote:
> Hi-
>
> I am trying to find a way to diffuse GIS data on a European map. I have a
> dataset consisting of particular locations scattered across Europe,
> along with magnitude and value information. I can plot these as discrete
> points with something like the following:
>
> "geocode" is a dataframe with four columns: LAT; LONG; MAGNITUDE;VALUE.
>
> library(maps)
> library(mapdata)
> map("worldHires", regions=c("Germany", "Belgium", "Netherlands"))
> points(geocode$LONG, geocode$LAT, cex=geocode$MAGNITUDE / 2500,
> col=rainbow(length(geocode$VALUE), start=0, end=.4)[rank(geocode$VALUE)])
>
> This gives me a map of Europe with my datapoints highlighted in two ways:
> magnitude is represented by the size of the point, and value is
> represented by the color.
>
> However, what I would really like is for there to be some sort of
> diffusion, such that instead of discrete points, the European map is
> covered in color so I can see more clearly whether there are regional
> patterns (something that will presumably look like this contour chart:
>   only
> on the European map).
>
> I have absolutely no idea where to start because I can't find a function
> that will allow me to diffuse the datapoints on a map.
>
> thank you for any help
> ldb
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>


-- 
Felix Andrews / 安福立
PhD candidate
Integrated Catchment Assessment and Management Centre
The Fenner School of Environment and Society
The Australian National University (Building 48A), ACT 0200
Beijing Bag, Locked Bag 40, Kingston ACT 2604
http://www.neurofractal.org/felix/
voice:+86_1051404394 (in China)
mobile:+86_13522529265 (in China)
mobile:+61_410400963 (in Australia)
xmpp:[EMAIL PROTECTED]
3358 543D AAC6 22C2 D336  80D9 360B 72DD 3E4C F5D8

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Linear Regression with slope equals 0

[E . N . D . Grew]

 Hi there, am trying to run a linear regression with a slope of 0.

 I have a dataset as follows

 t d
 1 303
 2 302
 3 304
 4 306
 5 307
 6 303

 I would like to test the significance that these points would lie on a
horizontal straight line.

 The standard regression lm(d~t) doesn't seem to allow the slope to be set.

 Any help very welcome.

 ed

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] pakages

[Uwe Ligges]

Safaa Najla wrote:
> Good Morning
> i am in doctorate.
> i want to ask for the packages of R,
> i tried to install it, 

What? R or packages?
How did you try it, please give us more details as suggested in the 
posting guide!

 > but it appears " not found".

Complete error message, please.

> each time i have to 
> delete and install again the programme R.

Why?


> Are there another method for installing the package?

Which one are we talking about? Or do you mean R itself? In that case, 
what is your OS etc.? Again, please read the posting guide.

Uwe Ligges


> thank you
> 
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] BDS test - results unclear to me

[Adrian Trapletti]

Hi Nicolas

Most likely because the higher embedding dimension uses less data 
points, even when applied to the lower dimension. Therefore, the actual 
numbers should be similar, but do not exactly match. But it's a very 
long time ago I last touched this code. So pls read bdstest.c from 
tseries to be sure.

Best regards
Adrian



> Subject:
> [R] BDS test - results unclear to me
> From:
> Nicolas Navet <[EMAIL PROTECTED]>
> Date:
> Mon, 13 Aug 2007 22:14:40 +0200
>
> To:
> r-help@stat.math.ethz.ch
>
>
> Hello,
>
> I would like to use the BDS test from the tseries package, but there 
> is something I don't understand in the results of the test. Let's say, 
> I want the BDS values for an embedding dimension equal to 2 :
>
> > bds.test(c, m = 2, eps = seq(0.5 * sd(c), 2 * sd(c), length = 
> 4),trace=FALSE);
> Here are the outputs:
> data:  c
> Embedding dimension =  2
> Epsilon for close points =  0.0097 0.0194 0.0291 0.0388
> Standard Normal =
>  [ 0.0097 ] [ 0.0194 ] [ 0.0291 ] [ 0.0388 ]
> [ 2 ]14.600613.900312.574511.4012
>
>
> Now for dimension 3 (m=3), we obtain:
>
> Standard Normal =
>  [ 0.0097 ] [ 0.0194 ] [ 0.0291 ] [ 0.0388 ]
> [ 2 ]14.554413.875812.554011.3897
> [ 3 ]20.914918.764016.156214.2518
>
> what I don't understand is why the values for embedding dimension 2 
> are not equal when BDS is computed with parameter m=2 and m=3, could 
> someone please explain that to me ? In the documentation, it is said 
> that "m is an integer indicating that the BDS test statistic is 
> computed for embedding dimensions 2, ..., m.", so why don't we get the 
> same result in both cases?
>
> Thank you for your help,
> Best regards,
>
> Nicolas
>

-- 
Adrian Trapletti
Wildsbergstrasse 31
8610 Uster
Switzerland

Phone :   +41 (0) 44 9945630
Mobile :  +41 (0) 76 3705631

Email :   [EMAIL PROTECTED]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] pakages

[gyadav]

please give details of the system, OS, etc which version of R you are 
using, how you are installing the packages

etc etc


---
  Regards,
Gaurav Yadav (mobile: +919821286118)
Assistant Manager, CCIL, Mumbai (India)
mailto:[EMAIL PROTECTED]
mailto:[EMAIL PROTECTED]
Profile: http://www.linkedin.com/in/gydec25
  Keep in touch and keep mailing :-)
slow or fast, little or too much




Safaa Najla <[EMAIL PROTECTED]> 
Sent by: [EMAIL PROTECTED]
08/14/2007 04:29 PM

To
R-help@stat.math.ethz.ch
cc

Subject
[R] pakages






Good Morning
i am in doctorate.
i want to ask for the packages of R,
i tried to install it, but it appears " not found". each time i have to 
delete and install again the programme R.
Are there another method for installing the package?
thank you

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide 
http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.




DISCLAIMER AND CONFIDENTIALITY CAUTION:\ \ This message and ...{{dropped}}

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] pakages

[Safaa Najla]

Good Morning
i am in doctorate.
i want to ask for the packages of R,
i tried to install it, but it appears " not found". each time i have to 
delete and install again the programme R.
Are there another method for installing the package?
thank you

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] GML with tweedie: AIC=NA

[Catarina Miranda]

 Thanks for your reply, but I still couldn't solve the problem.

I am using the package statmod, and I need the AIC because I want to use the
step function (I am modelling many species, so I would prefer to do the
step automatically).
I can't find the tweedie  package in the R packages list and I don't know
how to download it from the package source available in
*http://cran.r-project.org/src/contrib/Descriptions/tweedie.html
* .
However, after running the  tweedie.R file from that package source I am
able to use the dtweedie() function, but still I didn't figure out a way to
do the step (or to get the AIC from the glm command).

Thank you again for your help;

Catarina


On 13/08/07, Gordon Smyth <[EMAIL PROTECTED]> wrote:
>
> Dear Catarina,
>
> I prefer to leave the AIC value as NA for the tweedie GLM family
> because it takes extra time to compute and is only occasionally
> wanted. It's easy to compute the AIC yourself using the dtweedie()
> function of the tweedie package.
>
> Best wishes
> Gordon
>
> At 03:05 AM 14/08/2007, Catarina Miranda wrote:
> >Dear Gordon;
> >
> >I have also sent this email to R help mailing list, so I apologize
> >for duplicated mailing.
> >I am modelling densities of some species of birds, and I have a
> >problem with a great amount of zeros.
> >I have decided to try GLMs with the tweedie family, but in all the
> >models I have tried  I got an NA for the AIC value.
> >Just  to check the problem I've compared the a glm using the
> >Gaussian family with the identity link and a glm using the tweedie
> >family with var.power=0 and link.power=1. These are equal, as
> >expected, except the fact that the tweedie output gives me an NA for the
> AIC.
> >Could you help me with this problem?
> >Below you can find the two outputs I refer.
> >
> >Best Wishes;
> >
> >Catarina
> >
> > > summary(glm(formula=ACIN~DIST_REF+DIST_H2O+DIST_OST+
> > COTA+H2O_SUP+vasa,family=gaussian(link="identity")))
> >Call:glm(formula = ACIN ~ DIST_REF + DIST_H2O + DIST_OST + COTA
> >+ H2O_SUP + vasa, family = gaussian(link = "identity"))
> >Deviance
> >Residuals:   Min 1Q Median 3QMax
> >-0.112792   -0.042860  -0.021113  -0.006311   1.551824
> >Coefficients:  Estimate Std. Error t value
> >Pr(>|t|)  (Intercept)
> >-6.625e-02  5.454e-02  -1.215   0.2256  DIST_REF 3.581e-06
> >1.336e-05   0.268   0.7889  DIST_H2O-
> > 3.168e-05  1.527e-05  -2.074   0.0391
> >*DIST_OST-1.799e-05  1.953e-05  -0.921   0.3579  COTA
> >5.648e-04  2.470e-04   2.287   0.0230
> >*H2O_SUP -2.172e-04  3.994e-04  -0.544   0.5870  vasa
> >3.695e-02   4.573e-020.808   0.4199  ---Signif. codes:  0 '***'
> >0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
> >(Dispersion parameter for gaussian family taken to be 0.02151985)
> > Null deviance: 5.6028   on 257  degrees of freedomResidual
> > deviance: 5.4015  on 251  degrees of freedomAIC: -249.33
> >Number of Fisher Scoring iterations: 2
> >
> >
> > > summary(glm(formula=ACIN~DIST_REF+DIST_H2O+DIST_OST+
> > COTA+H2O_SUP+vasa,control=
> > glm.control(maxit=750),family=tweedie(var.power=0, link.power=1)))
> >Call:glm(formula = ACIN ~ DIST_REF + DIST_H2O + DIST_OST + COTA
> >+ H2O_SUP + vasa, family = tweedie( var.power = 0, link.power =
> >1), control = glm.control (maxit = 750))
> >Deviance
> >Residuals:   Min 1Q Median 3QMax
> >-0.112792  -0.042860  -0.021113  -0.006311   1.551824
> >Coefficients:  Estimate Std. Error t value
> >Pr(>|t|)  (Intercept) -
> >6.625e-02  5.454e-02  -1.215   0.2256  DIST_REF 3.581e-06
> >1.336e-05   0.268   0.7889  DIST_H2O-3.168e-05   1.527e-05
> >-2.074   0.0391
> >*DIST_OST-1.799e-05  1.953e-05  -0.921   0.3579  COTA
> >5.648e-04  2.470e-042.287   0.0230
> >*H2O_SUP -2.172e-04  3.994e-04  -0.544   0.5870  vasa
> >3.695e-02   4.573e-02   0.808   0.4199  ---Signif. codes:  0 '***'
> >0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
> >(Dispersion parameter for Tweedie family taken to be 0.02151985)
> > Null deviance: 5.6028  on 257  degrees of freedomResidual
> > deviance: 5.4015  on 251  degrees of freedomAIC: NA
> >Number of Fisher Scoring iterations: 2
> >
> >
>
>

[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Slack variable in OR

[Amir Safari]

 
   Hi 
  Thank you very much.  
   
  So sorry for the badly writing an example.
   
  Let assume we are modelling in a quadratic programming framework. So we have
   
  min 1/2xHx + fx subject to: Ax= 
Sent by: [EMAIL PROTECTED]   08/14/2007 02:37 PM 
To
  R-help@stat.math.ethz.ch   cc
Subject
  [R] Slack variable in OR
  



 
  
Hi dear R users,
  
@@@ although i have not understood your problemos clearly whatever i can help, 
i will tell 

 Is it basically correct that a problem is ( linearly on nonlinearly ) modeled 
so that the slack variable is bounded by an upper bound ? 
  
 If so, how it can be handled and coded practically ?
  
 for example:
  
 x1+ x2 =< b   so >  x1 + x2 + s=b
  
 s=b- x1 - x2
  
 b- x1 - x2 =< upper value
@@@ i assume that x1 and x2 are >= 0, which is usually the case 
@@@ then in that case if both x1=x2=0 then x3=50  which is its meximum 
value   
 But algorithms can not calculate b- x1 - x2 , because the matrices are not 
compatible. Or how can it be well coded ?
@@@ where are the matrices i am just asking :-)  
 Thank you so much.
 Amir

  
-
Boardwalk for $500? In 2007? Ha! 

[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.



DISCLAIMER AND CONFIDENTIALITY CAUTION:\ \ This message and ...{{dropped}}

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] graph dimensions default

[Simon Pickett]

Hi,

I would like to (if possible) set the default width and height for graphs
at the start of each session and have each new graphic device overwrite
the previous one.

I only know how to do this using windows(width=,height=...) which opens up
a new plotting device every time, so I end up with lots of graphs all over
the place until I get the one I want!

Thanks in advance,

Simon


Simon Pickett
PhD student
Centre For Ecology and Conservation
Tremough Campus
University of Exeter in Cornwall
TR109EZ
Tel 01326371852

http://www.uec.ac.uk/biology/research/phd-students/simon_pickett.shtml

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Help with npmc

[Paul Fenner]

I cant't seem to get npmc to make a comparison to a control level

>summary(npmc(brain), type="BF", control=1)
$`Data-structure`
  group.index class.level nobs
c   1   c   30
l   2   l   30
r   3   r   30

$`Results of the multiple Behrens-Fisher-Test`
  cmpeffect  lower.cl  upper.cl p.value.1s p.value.2s
1 1-2 0.643 0.4610459 0.8256208 0.08595894 0.14750647
2 1-3 0.444 0.2576352 0.6312537 0.99636221 0.75376639
3 2-3 0.328 0.1602449 0.4964218 1. 0.04476692

What elementary error am I making.
Thanks,
Paul

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] diffusing GIS data in maps

[Vladimir Eremeev]

Hi. 
You will find some useful information in the 
  http://r-spatial.sourceforge.net/ 

Particularly, Fig. 07 in the Graph gallery.

Package spmaps can be used to extract desired boundaries from the mapdata
and convert them to the format suitable for sp and others.


Lawrence D. Brenninkmeyer wrote:
> 
> Hi-
> 
> I am trying to find a way to diffuse GIS data on a European map. I have a
> dataset consisting of particular locations scattered across Europe,
> along with magnitude and value information. I can plot these as discrete
> points with something like the following:
> 
> "geocode" is a dataframe with four columns: LAT; LONG; MAGNITUDE;VALUE.
> 
> library(maps)
> library(mapdata)
> map("worldHires", regions=c("Germany", "Belgium", "Netherlands"))
> points(geocode$LONG, geocode$LAT, cex=geocode$MAGNITUDE / 2500,
> col=rainbow(length(geocode$VALUE), start=0, end=.4)[rank(geocode$VALUE)])
> 
> This gives me a map of Europe with my datapoints highlighted in two ways:
> magnitude is represented by the size of the point, and value is
> represented by the color.
> 
> However, what I would really like is for there to be some sort of
> diffusion, such that instead of discrete points, the European map is
> covered in color so I can see more clearly whether there are regional
> patterns (something that will presumably look like this contour chart:
>   only
> on the European map).
> 
> 

-- 
View this message in context: 
http://www.nabble.com/diffusing-GIS-data-in-maps-tf4265350.html#a12141809
Sent from the R help mailing list archive at Nabble.com.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Import of Access data via RODBC changes column name ("NO" to "Expr1014") and the content of the column

[Prof Brian Ripley]

On Tue, 14 Aug 2007, Maciej Hoffman-Wecker wrote:

>
> Dear all,
>
> I have some problems with importing data from an Access data base via
> RODBC to R. The data base contains several tables, which all are
> imported consecutively. One table has a column with column name "NO". If
> I run the code attached on the bottom of the mail I get no complain, but
> the column name (name of the respective vector of the data.frame) is
> "Expr1014" instead of "NO". Additionally the original column (type
> "text") containes "0"s and missings, but the imported column contains
> "0"s only (type "int"). If I change the column name in the Access data
> base to "NOx", the import works fine with the right name and the same
> data.
>
> Previously I generated a tiny Access data base which reproduced the
> problem. To be on the safe site I installed the latest version (2.5.1)
> and now the example works fine, but within my production process the
> error still remaines. An import into excel via ODBC works fine.
>
> So there is no way to figure it out whether this is a bug or a
> feature.-)

It's most likely an ODBC issue, but you have not provided a reproducible 
example.

> The second problem I have is that when I rerun "rm(list = ls(all = T));
> gc()" and the import several times I get the following error:
>
> Error in odbcTables(channel) : Calloc could not allocate (263168 of 1)
> memory
> In addition: Warning messages:
> 1: Reached total allocation of 447Mb: see help(memory.size) in:
> odbcQuery(channel, query, rows_at_time)
> 2: Reached total allocation of 447Mb: see help(memory.size) in:
> odbcQuery(channel, query, rows_at_time)
> 3: Reached total allocation of 447Mb: see help(memory.size) in:
> odbcTables(channel)
> 4: Reached total allocation of 447Mb: see help(memory.size) in:
> odbcTables(channel)
>
> which is surprising to me, as the first two statements should delete all

How do you _know _what they 'should' do?  That only deletes all objects in 
the workspace, not all objects in R, and not all memory blocks used by R.

Please do read ?"Memory-limits" for the possible reasons.

Where did '447Mb' come from?  If this machine has less than 2Gb of RAM, 
buy some more.


> objects and recover the memory. Is this only a matter of memory? Is
> there any logging that reduces the memory? Or is this issue connected to
> the upper problem?
>
> I added the code on the bottom - maybe there is some kind of misuse I
> lost sight of. Any hints are appreciated.
>
> Kind regards,
> Maciej
>
>> version
>   _
> platform   i386-pc-mingw32
> arch   i386
> os mingw32
> system i386, mingw32
> status
> major  2
> minor  5.1
> year   2007
> month  06
> day27
> svn rev42083
> language   R
> version.string R version 2.5.1 (2007-06-27)
>
>
> ## code
>
> get.table <- function(name, db, drop = NULL){
>  .con <- try(odbcConnectAccess(db), silent = T)
>  if(!inherits(.con, "RODBC")) return(.con)
>  ## exclude memo columns
>  .t <- try(sqlColumns(.con, name))
>  if(inherits(.t, "try-error")){close(.con); return(.t)}
>  .t <- .t[.t$"COLUMN_SIZE" < 255, "COLUMN_NAME"]
>  .t <- paste(.t, collapse = ",")
>  ## get table
>  .t <- paste("select", .t, "from", name)
>  .d <- try(sqlQuery(.con, .t), silent = T)
>  if(inherits(.d, "try-error")){close(.con); return(.d)}
>  .con <- try(close(.con), silent = T)
>  if(inherits(.con, "try-error")) return(.con)
>  .d <- .d[!names(.d) %in% drop]
>  return(.d)
> }
>
> get.alltables <- function(db){
>  .con <- try(odbcConnectAccess(db), silent = T)
>  if(!inherits(.con, "RODBC")) return(.con)
>  .tbls <- try(sqlTables(.con)[["TABLE_NAME"]])
>  if(inherits(.tbls, "try-error")){close(.con); return(.tbls)}
>  .con <- try(close(.con), silent = T)
>  if(inherits(.con, "try-error")) return(.con)
>  .tbls <- .tbls[-grep("^MSys", .tbls)]
>  .d <- lapply(seq(along = .tbls), function(.i){
>.d <-
>  try(get.table(.tbls[.i], db = db))
>return(invisible(.d))
>  })
>  names(.d) <- .tbls
>  .ok <- !sapply(.d, inherits, "try-error")
>  return(list(notdone = .d[!.ok], data = .d[.ok]))
> }
>
> library(RODBC)
>
> alldata <- get.alltables(db = "./myaccessdb.MDB")
>
> ## code end
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://

Re: [R] Slack variable in OR

[gyadav]

hi Amir

please see inline answers



Amir Safari <[EMAIL PROTECTED]> 
Sent by: [EMAIL PROTECTED]
08/14/2007 02:37 PM

To
R-help@stat.math.ethz.ch
cc

Subject
[R] Slack variable in OR






 
 
 Hi dear R users,
 
@@@ although i have not understood your problemos clearly whatever i can 
help, i will tell

  Is it basically correct that a problem is ( linearly on nonlinearly ) 
modeled so that the slack variable is bounded by an upper bound ? 
 
  If so, how it can be handled and coded practically ?
 
  for example:
 
  x1+ x2 =< b   so >  x1 + x2 + s=b
 
  s=b- x1 - x2
 
  b- x1 - x2 =< upper value
@@@ i assume that x1 and x2 are >= 0, which is usually the case
@@@ then in that case if both x1=x2=0 then x3=50  which is its meximum 
value 
  But algorithms can not calculate b- x1 - x2 , because the matrices are 
not compatible. Or how can it be well coded ?
 @@@ where are the matrices i am just asking :-) 
  Thank you so much.
  Amir

 
-
Boardwalk for $500? In 2007? Ha! 

 [[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide 
http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.




DISCLAIMER AND CONFIDENTIALITY CAUTION:\ \ This message and ...{{dropped}}

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Slack variable in OR

[Amir Safari]

 
   
 Hi dear R users,
   
  Is it basically correct that a problem is ( linearly on nonlinearly ) modeled 
so that the slack variable is bounded by an upper bound ?
   
  If so, how it can be handled and coded practically ?
   
  for example:
   
  x1+ x2 =< b   so >  x1 + x2 + s=b
   
  s=b- x1 - x2
   
  b- x1 - x2 =< upper value
   
  But algorithms can not calculate b- x1 - x2 , because the matrices are not 
compatible. Or how can it be well coded ?
   
  Thank you so much.
  Amir

   
-
Boardwalk for $500? In 2007? Ha! 

[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] vertically oriented color key in heatmaps

[Jim Lemon]

Rajarshi Guha wrote:
> Hi, I have some data which I was plotting using image(). I wanted to  
> add a vertical color key to the plot and I found that heatmap.2 in  
> gplots does let me add a color key. However, I was thinking of a  
> vertical bar with the color range rather than  the style that gplots  
> provides.
> 
> Is there any package (or code snippet) that would let me add a  
> vertical color key to an image() or heatmap plot?
> 
Hi Rajarshi,
Have a look at color.legend in the plotrix package.

Jim

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] RWeka cross-validation and Weka_control Parametrization

[Kurt Hornik]

> On Wed, 01 Aug 2007 10:52:02 +0200, Bjoern wrote:

> Hello,

>  I have two questions concerning the RWeka package:

>  1.) First question:
>  How can one perform a cross validation, -say 10fold- for a given 
> data set and given model ?

>  2.) Second question
>  What is the correct syntax for the parametrization of e.g. Kernel 
> classifiers interface
>m1 <- SMO(Species ~ ., data = iris, control = 
>  
> Weka_control(K="weka.classifiers.functions.supportVector.RBFKernel",G=0.1))
>m2 <- SMO(Species ~ ., data = iris, control = 
>  
> Weka_control(K="weka.classifiers.functions.supportVector.RBFKernel",G=1.0))

>> m1
>  SMO

>  Kernel used:
>  RBF kernel: K(x,y) = e^-(0.01* ^2)

>  ## should be: RBF kernel: K(x,y) = e^-(0.1* ^2)

> etc.

The answer for question 2 is surprisingly simple, but nevertheless took
me about half an hour to find:

  m2 <- SMO(Species ~ ., data = iris,
  control = Weka_control(K = 
"weka.classifiers.functions.supportVector.RBFKernel -G 2"))

gives

R> m2
SMO

Kernel used:
  RBF kernel: K(x,y) = e^-(2.0* ^2)

[Using Weka_control(K = ..., G = ...) passes the G option to SMO but not
RBFKernel.  The docs for SMO() say

 -K 
  The Kernel to use.
  (default: weka.classifiers.functions.supportVector.PolyKernel)

and one needs to remember Weka's command line style interface to realize
that this deparses into putting everything into a string for the K
option.]

This is of course not quite what R users would expect, and we'll try to
improve the Weka control mechanism so that specifying (Weka class)
options which require additional parameters becomes more convenient.

Best
-k

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] cov.unscaled in gls object

[Prof Brian Ripley]

This is what the vcov() generic is for.  You are asking for internal 
details from a different class ("summary.lm").

On Tue, 14 Aug 2007, Sven Garbade wrote:

> Hi list,
>
> can I extract the cov.unscaled ("the unscaled covariance matrix") from a
> gls fit (package nlme), like with summary.lm? Background: In a fixed
> effect meta analysis regression the standard errors of the coefficients
> can be computed as sqrt(diag(cov.unscaled)) where cov.unscaled is
> (X'WX). I try do do this with a gls-fit.

I don't think so: the 'unscaled' is a clue.  The vcov method is

> stats:::vcov.lm
function (object, ...)
{
 so <- summary.lm(object, corr = FALSE)
 so$sigma^2 * so$cov.unscaled
}

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Import of Access data via RODBC changes column name ("NO" to "Expr1014") and the content of the column

[Maciej Hoffman-Wecker]

Dear all,
 
I have some problems with importing data from an Access data base via
RODBC to R. The data base contains several tables, which all are
imported consecutively. One table has a column with column name "NO". If
I run the code attached on the bottom of the mail I get no complain, but
the column name (name of the respective vector of the data.frame) is
"Expr1014" instead of "NO". Additionally the original column (type
"text") containes "0"s and missings, but the imported column contains
"0"s only (type "int"). If I change the column name in the Access data
base to "NOx", the import works fine with the right name and the same
data. 
 
Previously I generated a tiny Access data base which reproduced the
problem. To be on the safe site I installed the latest version (2.5.1)
and now the example works fine, but within my production process the
error still remaines. An import into excel via ODBC works fine.

So there is no way to figure it out whether this is a bug or a
feature.-)
 
The second problem I have is that when I rerun "rm(list = ls(all = T));
gc()" and the import several times I get the following error:
 
Error in odbcTables(channel) : Calloc could not allocate (263168 of 1)
memory
In addition: Warning messages:
1: Reached total allocation of 447Mb: see help(memory.size) in:
odbcQuery(channel, query, rows_at_time) 
2: Reached total allocation of 447Mb: see help(memory.size) in:
odbcQuery(channel, query, rows_at_time) 
3: Reached total allocation of 447Mb: see help(memory.size) in:
odbcTables(channel) 
4: Reached total allocation of 447Mb: see help(memory.size) in:
odbcTables(channel) 

which is surprising to me, as the first two statements should delete all
objects and recover the memory. Is this only a matter of memory? Is
there any logging that reduces the memory? Or is this issue connected to
the upper problem?
 
I added the code on the bottom - maybe there is some kind of misuse I
lost sight of. Any hints are appreciated.
 
Kind regards,
Maciej
 
> version
   _   
platform   i386-pc-mingw32 
arch   i386
os mingw32 
system i386, mingw32   
status 
major  2   
minor  5.1 
year   2007
month  06  
day27  
svn rev42083   
language   R   
version.string R version 2.5.1 (2007-06-27)

 
## code
 
get.table <- function(name, db, drop = NULL){
  .con <- try(odbcConnectAccess(db), silent = T)
  if(!inherits(.con, "RODBC")) return(.con)
  ## exclude memo columns
  .t <- try(sqlColumns(.con, name))
  if(inherits(.t, "try-error")){close(.con); return(.t)}
  .t <- .t[.t$"COLUMN_SIZE" < 255, "COLUMN_NAME"]
  .t <- paste(.t, collapse = ",")
  ## get table
  .t <- paste("select", .t, "from", name)
  .d <- try(sqlQuery(.con, .t), silent = T)
  if(inherits(.d, "try-error")){close(.con); return(.d)}
  .con <- try(close(.con), silent = T)
  if(inherits(.con, "try-error")) return(.con)
  .d <- .d[!names(.d) %in% drop]
  return(.d)
}
 
get.alltables <- function(db){
  .con <- try(odbcConnectAccess(db), silent = T)
  if(!inherits(.con, "RODBC")) return(.con)
  .tbls <- try(sqlTables(.con)[["TABLE_NAME"]])
  if(inherits(.tbls, "try-error")){close(.con); return(.tbls)}
  .con <- try(close(.con), silent = T)
  if(inherits(.con, "try-error")) return(.con)
  .tbls <- .tbls[-grep("^MSys", .tbls)]
  .d <- lapply(seq(along = .tbls), function(.i){
.d <-
  try(get.table(.tbls[.i], db = db)) 
return(invisible(.d))
  })
  names(.d) <- .tbls
  .ok <- !sapply(.d, inherits, "try-error")
  return(list(notdone = .d[!.ok], data = .d[.ok]))
}
 
library(RODBC)
 
alldata <- get.alltables(db = "./myaccessdb.MDB")
 
## code end

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] R 2.5.1 configure problem

[Peter Dalgaard]

Andreas Hey wrote:
> Hi,
>
>  
>
> I have follwoing problem:
>
>  
>
> I will install R-2.5.1 on a Linux Maschine (64Bit) 

Which? (CPU and OS, please. There are about four likely possibilities,
half a dozen less likely ones...)

> and I will use the R-GUI
> JGR (Jaguar)
>
>  
>
> SO I make following steps:
>
>  
>
> ./configure --with-gnu-ld --enable-R-shlib VAR=fPIC VAR=TCLTK_LIBS
>
>   
What are those options supposed to be good for??? You appear to be
setting VAR twice, and I don't recall VAR as anything used by configure.
And it is unlikely that a Linux system would use anything but GNU ld by
default. Did configure terminate succesfully??? What did  the output
summary say?

> make 
>
>  
>
> I become following error messages:
>
>   
(that's not how to translate "Ich bekomme"...)
>  
>
> .
>   
Something must have gone before this! A linker error perhaps?
> Entering directory `/root/R-2.5.1/tests'
>
> make[2]: Entering directory `/root/R-2.5.1/tests'
>
> make[3]: Entering directory `/root/R-2.5.1/tests/Examples'
>
> make[4]: Entering directory `/root/R-2.5.1/tests/Examples'
>
> make[4]: `Makedeps' is up to date.
>
> make[4]: Leaving directory `/root/R-2.5.1/tests/Examples'
>
> make[4]: Entering directory `/root/R-2.5.1/tests/Examples'
>
> make[4]: *** No rule to make target `../../lib/libR.so', needed by
> `base-Ex.Rout'.  Stop.
>   
Did you really run "make". This looks like "make check" output.

> make[4]: Leaving directory `/root/R-2.5.1/tests/Examples'
>
> make[3]: *** [test-Examples-Base] Error 2
>
> make[3]: Leaving directory `/root/R-2.5.1/tests/Examples'
>
> make[2]: *** [test-Examples] Error 2
>
> make[2]: Leaving directory `/root/R-2.5.1/tests'
>
> make[1]: *** [test-all-basics] Error 1
>
> make[1]: Leaving directory `/root/R-2.5.1/tests'
>
> make: *** [check] Error 2
>
>  
>
> Make check all - I become following messages:
>
>  
>
> Entering directory `/root/R-2.5.1/tests'
>
> make[2]: Entering directory `/root/R-2.5.1/tests'
>
> make[3]: Entering directory `/root/R-2.5.1/tests/Examples'
>
> make[4]: Entering directory `/root/R-2.5.1/tests/Examples'
>
> make[4]: `Makedeps' is up to date.
>
> make[4]: Leaving directory `/root/R-2.5.1/tests/Examples'
>
> make[4]: Entering directory `/root/R-2.5.1/tests/Examples'
>
> make[4]: *** No rule to make target `../../lib/libR.so', needed by
> `base-Ex.Rout'.  Stop.
>
> make[4]: Leaving directory `/root/R-2.5.1/tests/Examples'
>
> make[3]: *** [test-Examples-Base] Error 2
>
> make[3]: Leaving directory `/root/R-2.5.1/tests/Examples'
>
> make[2]: *** [test-Examples] Error 2
>
> make[2]: Leaving directory `/root/R-2.5.1/tests'
>
> make[1]: *** [test-all-basics] Error 1
>
> make[1]: Leaving directory `/root/R-2.5.1/tests'
>
> make: *** [check] Error 2
>
>  
>
>  
>
> Can you help me?
>
>  
>
> With best regards 
>
>  
>
> Andreas Hey
>
>  
>
> _
>
> Tel:030/2093-1463
>
> Email:   [EMAIL PROTECTED]
>
>  
>
>
>   [[alternative HTML version deleted]]
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>   


-- 
   O__   Peter Dalgaard Øster Farimagsgade 5, Entr.B
  c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K
 (*) \(*) -- University of Copenhagen   Denmark  Ph:  (+45) 35327918
~~ - ([EMAIL PROTECTED])  FAX: (+45) 35327907

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] binomial simulation

[sigalit mangut-leiba]

Thank you,
I'm trying to run the joint probabilty:

C(N,k)*p^k*(1-p)^(N-k)*C(k,m)*q^m*(1-q)^(k-m)

and get the error: Error in C(N, k) : object not interpretable as a factor

so I tried the long way:

gamma(N+1)/(gamma(k+1)*(gamma(N-k)))

and the same with k, and got the error:

1: value out of range in 'gammafn' in: gamma(N + 1)
2: value out of range in 'gammafn' in: gamma(N - k) 

Do you know why it's not working?

Thanks again,

Sigalit.



On 8/14/07, Moshe Olshansky <[EMAIL PROTECTED]> wrote:
>
> As I understand this,
> P(T+ | D-)=1-P(T+ | D+)=0.05
> is the probability not to detect desease for a person
> at ICU who has the desease. Correct?
>
> What I asked was whether it is possible to mistakenly
> detect the desease for a person who does not have it?
>
> Assuming that this is impossible the formula is below:
>
> If there are N patients, each has a probability p to
> have the desease (p=0.6 in your case) and q is the
> probability to detect the desease for a person who has
> it (q = 0.95 for ICU and q = 0.8 for a regular unit),
> then
>
> P(k have the desease AND m are detected) =
> P(k have the desease)*P(m are detected / k have the
> desease) =
> C(N,k)*p^k*(1-p)^(N-k)*C(k,m)*q^m*(1-q)^(k-m)
> where C(a,b) is the Binomial coefficient "a above b" -
> the number of ways to choose b items out of a (when
> the order does not matter). You of course must assume
> that N >= k >= m >= 0 (otherwise the probability is
> 0).
>
> To generate such pairs (k infected and m detected) you
> can do the following:
>
> k <- rbinom(N,1,p)
> m <- rbinom(k,1,q)
>
> Regards,
>
> Moshe.
>
> --- sigalit mangut-leiba <[EMAIL PROTECTED]> wrote:
>
> > Hi,
> > The probability of false detection is: P(T+ |
> > D-)=1-P(T+ | D+)=0.05.
> > and I want to find the joint probability
> > P(T+,D+)=P(T+|D+)*P(D+)
> > Thank you for your reply,
> > Sigalit.
> >
> >
> > On 8/13/07, Moshe Olshansky <[EMAIL PROTECTED]>
> > wrote:
> > >
> > > Hi Sigalit,
> > >
> > > Do you want to find the probability P(T+ = t AND
> > D+ =
> > > d) for all the combinations of t and d (for ICU
> > and
> > > Reg.)?
> > > Is the probability of false detection (when there
> > is
> > > no disease) always 0?
> > >
> > > Regards,
> > >
> > > Moshe.
> > >
> > > --- sigalit mangut-leiba <[EMAIL PROTECTED]>
> > wrote:
> > >
> > > > hello,
> > > > I asked about this simulation a few days ago,
> > but
> > > > still i can't get what i
> > > > need.
> > > > I have 2 units: icu and regular. from icu I want
> > to
> > > > take 200 observations
> > > > from binomial distribution, when probability for
> > > > disease is: p=0.6.
> > > > from regular I want to take 300 observation with
> > the
> > > > same probability: p=0.6
> > > > .
> > > > the distribution to detect disease when disease
> > > > occurred- *for someone from
> > > > icu* - is: p(T+ | D+)=0.95.
> > > > the distribution to detect disease when disease
> > > > occurred- *for someone from
> > > > reg.unit* - is: p(T+ | D+)=0.8.
> > > > I want to compute the joint distribution for
> > each
> > > > unit: p(T+,D+) for icu,
> > > > and the same for reg.
> > > > I tried:
> > > >
> > > > pdeti <- 0
> > > >
> > > > pdetr <- 0
> > > >
> > > > picu <- pdeti*.6
> > > >
> > > > preg <- pdetr*.6
> > > >
> > > > dept <- c("icu","reg")
> > > >
> > > > icu <- rbinom(200, 1, .6)
> > > >
> > > > reg <- rbinom(300, 1, .6)
> > > >
> > > > for(i in 1:300) {
> > > >
> > > > if(dept=="icu") pdeti==0.95
> > > >
> > > > if (dept=="reg") pdetr== 0.80
> > > >
> > > > }
> > > >
> > > > print(picu)
> > > >
> > > > print(preg)
> > > >
> > > > and got 50 warnings:
> > > >
> > > > the condition has length > 1 and only the first
> > > > element will be used in: if
> > > > (dept == "icu") pdeti == 0.95
> > > > the condition has length > 1 and only the first
> > > > element will be used in: if
> > > > (dept == "reg") pdetr == 0.8
> > > >
> > > > I would appreciate any suggestions,
> > > >
> > > > thank you,
> > > >
> > > > Sigalit.
> > > >
> > > >   [[alternative HTML version deleted]]
> > > >
> > > > __
> > > > R-help@stat.math.ethz.ch mailing list
> > > > https://stat.ethz.ch/mailman/listinfo/r-help
> > > > PLEASE do read the posting guide
> > > > http://www.R-project.org/posting-guide.html
> > > > and provide commented, minimal, self-contained,
> > > > reproducible code.
> > > >
> > >
> > >
> >
> >   [[alternative HTML version deleted]]
> >
> > __
> > R-help@stat.math.ethz.ch mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide
> > http://www.R-project.org/posting-guide.html
> > and provide commented, minimal, self-contained,
> > reproducible code.
> >
>
>

[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/l

[R] Mixed effects and hierarchical error terms

[Richard Gunton]

Please could someone help me understand the difference between the following two
models?

lme (Y ~ X, random= ~1|G1/G2)
aov (Y ~ X + Error (G1/G2))

Since random effects can be considered as additional error terms, I would have
thought these two model formulations should amount to the same analysis, but
the outputs (using suitable data, with X, G1 and G2 being factors) look very
different.

I note that the lme() model can be fitted even where the levels of X do not vary
within levels of G2, whereas the aov() model output is accompanied by a warning
message ("Error() model is singular...") in this case.

Thanks,

Richard.


--
Richard Gunton
PhD student - Ecology and Evolution group
School of Biology, University of Leeds, LS2 9JT, UK
Room 10.16, Miall Building   Tel: 0113 3432825

http://www.personal.leeds.ac.uk/~bgyrg

~ Opinions expressed in this message are not attributable to the University of
Leeds ~

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] R 2.5.1 configure problem

[Andreas Hey]

Hi,

 

I have follwoing problem:

 

I will install R-2.5.1 on a Linux Maschine (64Bit) and I will use the R-GUI
JGR (Jaguar)

 

SO I make following steps:

 

./configure --with-gnu-ld --enable-R-shlib VAR=fPIC VAR=TCLTK_LIBS

make 

 

I become following error messages:

 

.

Entering directory `/root/R-2.5.1/tests'

make[2]: Entering directory `/root/R-2.5.1/tests'

make[3]: Entering directory `/root/R-2.5.1/tests/Examples'

make[4]: Entering directory `/root/R-2.5.1/tests/Examples'

make[4]: `Makedeps' is up to date.

make[4]: Leaving directory `/root/R-2.5.1/tests/Examples'

make[4]: Entering directory `/root/R-2.5.1/tests/Examples'

make[4]: *** No rule to make target `../../lib/libR.so', needed by
`base-Ex.Rout'.  Stop.

make[4]: Leaving directory `/root/R-2.5.1/tests/Examples'

make[3]: *** [test-Examples-Base] Error 2

make[3]: Leaving directory `/root/R-2.5.1/tests/Examples'

make[2]: *** [test-Examples] Error 2

make[2]: Leaving directory `/root/R-2.5.1/tests'

make[1]: *** [test-all-basics] Error 1

make[1]: Leaving directory `/root/R-2.5.1/tests'

make: *** [check] Error 2

 

Make check all - I become following messages:

 

Entering directory `/root/R-2.5.1/tests'

make[2]: Entering directory `/root/R-2.5.1/tests'

make[3]: Entering directory `/root/R-2.5.1/tests/Examples'

make[4]: Entering directory `/root/R-2.5.1/tests/Examples'

make[4]: `Makedeps' is up to date.

make[4]: Leaving directory `/root/R-2.5.1/tests/Examples'

make[4]: Entering directory `/root/R-2.5.1/tests/Examples'

make[4]: *** No rule to make target `../../lib/libR.so', needed by
`base-Ex.Rout'.  Stop.

make[4]: Leaving directory `/root/R-2.5.1/tests/Examples'

make[3]: *** [test-Examples-Base] Error 2

make[3]: Leaving directory `/root/R-2.5.1/tests/Examples'

make[2]: *** [test-Examples] Error 2

make[2]: Leaving directory `/root/R-2.5.1/tests'

make[1]: *** [test-all-basics] Error 1

make[1]: Leaving directory `/root/R-2.5.1/tests'

make: *** [check] Error 2

 

 

Can you help me?

 

With best regards 

 

Andreas Hey

 

_

Tel:030/2093-1463

Email:   [EMAIL PROTECTED]

 


[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] makeSOCKcluster

[Michael Janis]

Hi,

 

I am attempting to implement a mixed (windows/linux) snow sockets
parallelism in R, but am running into difficulties similar to a post made
Aug 31, 2006 under the same subject heading.  I feel like I may be one or
two non-obvious steps away from getting it all working, but I'm stuck.  If
anyone can shed some light on this (I believe Prof. Tierney stated that he
has successfully run a snow master on Windows with slave nodes on Linux
using ssh.exe through Cygwin, which is exactly what I am attempting), I'd be
most grateful.

 

SNOW master on Windows information:

WinXP 32bit

R-2.5.1

I built a windows package (.zip) version of snow-0.2-3 for R using Rtools
and installed it without trouble.

Cygwin including ssh as an exe and in the path

 

Linux slave nodes information:

ROCKS compute nodes, each 64bit (I mention this difference between the win
and linux platforms out of desperation but I don't think it's an issue since
snow should be agnostic).  (btw, the SOCK version of snow here is just for
me to get used to parallelism in R before attempting MPI).

R-2.5.1, local directory install (no su, global R install is older).

Snow installed to R-2.5.1 instance

Unattended ssh authentication through public/private key pairs.  Each node
allows since they are all NFS.

 

Rgui commands:

>library(snow)

>cl<-makeCluster(c("ip.path.to.node1","ip.path.to.node2"),type="SOCK")

Result:

R becomes unresponsive and has to be forcibly closed.  Sometimes before that
happens message stating env: C:\pathToR\library\snow\RSOCKnode.sh no file
found appears (if I make R redraw the interface screen)

Troubleshooting step1:

Interspersed print commands within C:\localRinstall\library\snow\R\snow
file.

Result of troubleshooting step1:

It appears that this snow file hangs on line:
system(paste(rshcmd, "-l", user, machine, "env", env, script))

Resolution of troubleshooting step1:

  Manually attempted paste commands in R:

  >system("ssh -l me ip.path.to.node1 ls")

Result: gets ls of remote machine; so the script path must be
incorrect and the offending command

  Located RSOCKnode.sh on remote machine (for all remote machines, since
it's NFS)

Inserted line into C:\localRinstall\library\snow\R\snow file
just before the offending line system(paste(rshcmd, "-l", user, machine,
"env", env, script)).  That inserted line hard-codes the path of the script
and the script name like so:
script<-"/rootpath/home/cluster/me/R-2.5.1/library/snow/RSOCKnode.sh".  This
worked, and the next print statement following the system call in the snow
file now prints to the screen.  But again R hangs, this time at the line:
con <- socketConnection(port = port, server=TRUE, blocking=TRUE,
open="a+b").

Troubleshooting step2:

  Attempted manual socketConnection via R commands:

  > con<-socketConnection("ip.path.to.node1",port=22)

> showConnections()

description classmode text   isopen   can read can
write

3  "->ip.path.to.node1" "socket" "a+" "text" "opened" "yes""yes"

*** so the socketConnection does work, but in SNOW it hangs at this last
command, and I'm completely stuck.

 

As a parallel test, I have installed SNOW on the R-2.5.1 instance on the
linux cluster.  If I access the head node and launch R, I can use SNOW with
sockets successfully (from the linux head node as master to the linux
compute nodes as slaves). However, running R more or less continually on the
head node of a cluster is bad form.  Ideally I would like to run a windows
snow master to the linux node slaves, but I can't seem to get past this
point.

 

For reference, I include the original post below.  I am stuck at one command
in the snow file past where this last post had problems.

 

I hope this is clear.  I am intent on solving this problem, so feel free to
ask questions if you have feedback but my description is not clear.  I
really appreciate any help!

 

Yours,

 

Michael Janis

UCLA Bioinformatics

 

*transcript from original post entitles "makeSOCKcluster" follows*

 

 

makeSOCKcluster

Click to flag this post

 

by Hrishikesh Rajpathak Aug 31, 2006; 10:39pm :: Rate this Message: - Use
ratings to moderate (?)

 

Reply | Reply to Author | View Threaded | Show Only this Message

Hi,

 

  I am a newbie to R and trying to implement parallelism in R. I am
currently using R-2.3.1, and Cygwin to run R on Windows xp.

 

  ssh and all are working fine,

 

  When I try to create a socket connection as

 

   makeSOCKcluster(c("localhost","localhost")),

 

  it just waits for the other prcess on localhost to get created and
respond. But this other process is not created.

 

  To debug, I put print statements in the "snow " file in library\snow\r  of
R after every statement that comes under Socket Implementation. I  realized
that it does the execution till

 

  system(a<-paste(rshcmd, "-l", user, machine, "env", env, script))

 

  part and then it goes in wait stat

[R] LimmaGUI and Design Matrix

[tee` suesean]

Hi all,
 
I'm working on microarray and currently analyzing the microrarray data 
using limmaGUI. Loop design has been applied in this experiment. This is a 2X2 
factorial experiment where there are control and treatment at 2 different time 
points, week 6 and 9. The experimental design is almost the same as the 
limmaGUI work example: Weaver Data set. I would like to look at the effect of 
time, treatment and interaction. I have tried to follow the limmaGUI help 
example to construct the design matrix. The target file and design matrix are 
as follow: (C - control, M - treated)
 
SlideNumberFileNameCy3Cy5
M6.C6tss37-normnew.gprM6C6
C6.C9tss5-new.gprC6C9
C9.M9tss9-new.gprC9M9
M9.M6tss24-normnew.gprM9M6
C6.M6tss31-normnew.gprC6M6
M6.M9tss2-normnew.gprM6M9
M9.C9tss12-normnew.gprM9C9
C9.C6tss25-normnew.gprC9C6

 
SlideNumberTreatmentTimeInteraction
M6.C6-100
C6.C9010
C9.M9101
M9.M60-1-1
C6.M6100
M6.M9011
M9.C9-10-1
C9.C60-10

 
According to the design matrix above, may I know what is the meaning of time 
effect? Does it look at the differential expression of Control sample at 
different time? How about treatment effect? Is it refer to the differential 
gene expression for Treated sample at week 6? For the effect of interaction, is 
it refer to both time and treatment refer, which is refer to M9?
 
Besides that, I can not figure out the second example of design matrix, Mt for 
Weaver data set in LimmaGUI help file. Is the design matrix looking at the 
effect of different time point for mutant sample? Mt11 and Mt21 refer to the 
effect of mutant at day 11 and day 21 alone? What is the logically theory 
behind to asign 1 or -1 to slides P21WT.P11WT for effect of Mt11 and Mt21? Does 
it has anything to do with the vector? Or we will subtract it with Wt at both 
time point in order for us to look at the effect of Mt11 and Mt21?
 
I use Lucidea Universal ScoreCard as my spike in control. Do I need to flag 
them out before I carry out the analysis? Although I have include a Spot Type 
file, LimmaGUI still read all the spots as genes. 
 
Thank you 
 
 
Warm regards,
Tee Sue Sean
Department of Cell and Molecular Biology,
Faculty of Biotechnology and Biomolecular Sciences,
University Putra Malaysia, 
43400 Serdang, Selangor, Malaysia.


   


eChase.

[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] invert 160000x160000 matrix

[Prof Brian Ripley]

On Tue, 14 Aug 2007, Patnaik, Tirthankar  wrote:

> A variety of tricks would need to be used to invert a matrix of this 
> size. If there are any other properties of the matrix that you know 
> (symmetric, positive definite, etc, sparse) then they could be useful 
> too. You could partition the matrix first, then use an in-place inverse 
> technique for each block to individually calculate the blocks-inverses, 
> then combine to get the inverse of the initial matrix. Again, if the 
> implementation is actually solving an Ax-B = 0 system of equations, then 
> there are specific methods for these too, like an LU decomp, for 
> instance. You might also want to check out some texts for this, like the 
> Numerical Recipes.

> How's the matrix stored right now?

Well, not in R as a matrix: see ?"Memory-limits".  It is about 12x larger 
than the largest possible matrix in R.

>
> Best,
> -Tir
>
> Tirthankar Patnaik
> India Strategy
> Citigroup Investment Research
> +91-22-6631 9887
>
>> -Original Message-
>> From: [EMAIL PROTECTED]
>> [mailto:[EMAIL PROTECTED] On Behalf Of Moshe Olshansky
>> Sent: Tuesday, August 14, 2007 6:40 AM
>> To: Paul Gilbert; Jiao Yang
>> Cc: r-help@stat.math.ethz.ch
>> Subject: Re: [R] invert 16x16 matrix
>>
>> While inverting the matrix may be a problem, if you need to
>> solve an equation A*x = b you do not need to invert A, there
>> exist iterative methods which do need A or inv(A) - all you
>> need to provide is a function that computes A*x for an
>> arbitrary vector x.
>> For such a large matrix this may be slow but possible.
>>
>> --- Paul Gilbert <[EMAIL PROTECTED]>
>> wrote:
>>
>>> I don't think you can define a matrix this large in R, even if you
>>> have the memory. Then, of course, inverting it there may be other
>>> programs that have limitations.
>>>
>>> Paul
>>>
>>> Jiao Yang wrote:

 Can R invert a 16x16 matrix with all
>>> positive numbers?  Thanks a lot!

>

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] cov.unscaled in gls object

[Sven Garbade]

Hi list,

can I extract the cov.unscaled ("the unscaled covariance matrix") from a
gls fit (package nlme), like with summary.lm? Background: In a fixed
effect meta analysis regression the standard errors of the coefficients
can be computed as sqrt(diag(cov.unscaled)) where cov.unscaled is
(X'WX). I try do do this with a gls-fit. 

Thanks, Sven

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.