Re: [R] Problem Reading SPlus Dump Into R - Spaces Embedded in Data

2005-12-29 Thread allan miller


Peter Dalgaard wrote:

>allan miller <[EMAIL PROTECTED]> writes:
>
>  
>
>>Hello,
>>
>>I'm trying to source() an SPlus 6.x file created using dump(..., 
>>oldStyle=T) into R (version 2.01) as using the following instructions:
>>
>>
>>
>>>*If you have access to S-PLUS, it is usually more reliable to |dump| 
>>>the object(s) in S-PLUS and |source| the dumpfile in R. For S-PLUS 5.x 
>>>and 6.x you may need to use |dump(..., oldStyle=T)|, and to read in 
>>>very large objects it may be preferable to use the dumpfile as a batch 
>>>script rather than use the |source| function.*
>>>  
>>>
>>(from "R Data Import/Export," pg. 15)
>>
>>An example:
>>
>> > source("testdump")
>>Error in parse(file, n, text, prompt) : syntax error on line 1895
>>
>>where the data on line 1895 - and other lines causing this - have 
>>embedded spaces, such as the following:
>>
>>
>>[line 1895]  Johnson Partners LLC
>>
>>
>>I can't seem to find any options for either the SPlus dump, or R 
>>source(), that relate to this problem.  Any suggestions for how to 
>>either dump or source files containing data with embedded spaces?
>>
>>
>
>A bit more context might be helpful. What's in lines surrounding 1895?
>Can you show a simple S-PLUS object displaying the behaviour? What
>happens if you dput() the object? Will S-PLUS itself restore the file? 
>
>  
>
Unfortunately, I don't have access to S-PLUS :'( , the S-PLUS file dump 
was provided to me to load in R.  Here are the lines in the S-PLUS dump 
surrounding 1895:

> 1892 .Label
>1893 character
>1894 1
>1895 Johnson Partners LLC
>1896 class
>1897 character
>1898 1
>1899 factor
>1900 Protocol

The problem is with the embedded spaces (whitespace?) characters in 
1895.  If I remove the spaces, i.e., change it to:

JohnsonPartnersLLC

the line is successfully loaded, and the next error that comes up is 
another Label with embedded spaces.

Thanks for your help.

[[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


Re: [R] How to combine two lists?

2005-12-29 Thread Gabor Grothendieck
c(list(a=1,b=2), list(c=3,d=4))

or

append(list(a=1,b=2), list(c=3,d=4))


On 12/30/05, Xiao Shi <[EMAIL PROTECTED]> wrote:
>
> Hi every,
> I have two lists,one has 10 elements,another has 20 elements.And each
> element in these two lists has a unique name.I want to know  how to
> combine
> these two lists, so i can index the elements by their names.
>
> Thanks.
> Jiantao Shi
>
>[[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
>

[[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


[R] How to combine two lists?

2005-12-29 Thread Xiao Shi
Hi every,
I have two lists,one has 10 elements,another has 20 elements.And each
element in these two lists has a unique name.I want to know  how to combine
these two lists, so i can index the elements by their names.

Thanks.
Jiantao Shi

[[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


Re: [R] ESS and Emacs

2005-12-29 Thread paul sorenson
I tried it with XEmacs on Win XP and I had to install ESS separately.

Mark Leeds wrote:
> I have been using the document written by John Fox titled
> Sn Introduction to ESS + XEmacs for Windows
> Users of R.
>  
> It's a very nice document and
> I went through it carefully but I got
> an error when I finished it and launched XEmacs.
>  
> The error is "cannot open load file : ess-site".
>  
> So, I did more investigation
> and it seems like there is a folder
>  
> Program Files/XEmacs/xemacs-packages/etc/ess
>  
> that should have been created when
> I downloaded XEmacs but it doesn't exist. 
> There is a folder called efs ( rather than ess ) but
> I don't think that's it.
>  
> I tried installing XEmacs again but the
> same thing happened.
>  
> Does anyone know anything about this ?
> I am using Windows NT.
> Thanks.
>  
>  Mark
> 
> 
> **
> This email and any files transmitted with it are confidentia...{{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
>

__
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


[R] lme X lmer results

2005-12-29 Thread John Maindonald
Surely there is a correct denominator degrees of freedom if the design
is balanced, as Ronaldo's design seems to be. Assuming that he has
specified the design correctly to lme() and that lme() is getting the df
right, the difference is between 2 df and 878 df.  If the t-statistic  
for the
second level of Xvar had been 3.0 rather than 1.1, the difference
would be between a t-statistic equal to 0.095 and 1e-6.  In a design
where there are 10 observations on each experimental unit, and all
comparisons are at the level of experimental units or above, df for
all comparisons will be inflated by a factor of at least 9.

Rather than giving df that for the comparison(s) of interest may be
highly inflated, I'd prefer to give no degrees of freedom at all, & to
encourage users to work out df for themselves if at all possible.
If they are not able to do this, then mcmcsamp() is a good alternative,
and may be the way to go in any case.  This has the further advantage
of allowing assessments in cases where the relevant distribution is
hard to get at. I'd think a warning in order that the df are upper  
bounds,
and may be grossly inflated.

Incidentally, does mcmcsamp() do its calculations pretty well
independently of the lmer results?

John Maindonald.

On 29 Dec 2005, at 10:00 PM, [EMAIL PROTECTED] wrote:

> From: Douglas Bates <[EMAIL PROTECTED]>
> Date: 29 December 2005 5:59:07 AM
> To: "Ronaldo Reis-Jr." <[EMAIL PROTECTED]>
> Cc: R-Help 
> Subject: Re: [R] lme X lmer results
>
>
> On 12/26/05, Ronaldo Reis-Jr. <[EMAIL PROTECTED]> wrote:
>> Hi,
>>
>> this is not a new doubt, but is a doubt that I cant find a good  
>> response.
>>
>> Look this output:
>>
>>> m.lme <- lme(Yvar~Xvar,random=~1|Plot1/Plot2/Plot3)
>>
>>> anova(m.lme)
>> numDF denDF  F-value p-value
>> (Intercept) 1   860 210.2457  <.0001
>> Xvar1 2   1.2352  0.3821
>>> summary(m.lme)
>> Linear mixed-effects model fit by REML
>>  Data: NULL
>>   AIC  BIClogLik
>>   5416.59 5445.256 -2702.295
>>
>> Random effects:
>>  Formula: ~1 | Plot1
>> (Intercept)
>> StdDev: 0.000745924
>>
>>  Formula: ~1 | Plot2 %in% Plot1
>> (Intercept)
>> StdDev: 0.000158718
>>
>>  Formula: ~1 | Plot3 %in% Plot2 %in% Plot1
>> (Intercept) Residual
>> StdDev: 0.000196583 5.216954
>>
>> Fixed effects: Yvar ~ Xvar
>>Value Std.Error  DF  t-value p-value
>> (Intercept)2.3545454 0.2487091 860 9.467066  0.
>> XvarFactor20.3909091 0.3517278   2 1.111397  0.3821
>>
>> Number of Observations: 880
>> Number of Groups:
>>  Plot1   Plot2 %in% Plot1
>>  4  8
>>Plot3 %in% Plot2 %in% Plot1
>> 20
>>
>> This is the correct result, de correct denDF for Xvar.
>>
>> I make this using lmer.
>>
>>> m.lmer <- lmer(Yvar~Xvar+(1|Plot1)+(1|Plot1:Plot2)+(1|Plot3))
>>> anova(m.lmer)
>> Analysis of Variance Table
>>Df Sum Sq Mean Sq  Denom F value Pr(>F)
>> Xvar  1  33.62   33.62 878.00  1.2352 0.2667
>>> summary(m.lmer)
>> Linear mixed-effects model fit by REML
>> Formula: Yvar ~ Xvar + (1 | Plot1) + (1 | Plot1:Plot2) + (1 | Plot3)
>>  AIC BIClogLik MLdeviance REMLdeviance
>>  5416.59 5445.27 -2702.295   5402.698  5404.59
>> Random effects:
>>  GroupsNameVariance   Std.Dev.
>>  Plot3 (Intercept) 1.3608e-08 0.00011665
>>  Plot1:Plot2   (Intercept) 1.3608e-08 0.00011665
>>  Plot1 (Intercept) 1.3608e-08 0.00011665
>>  Residual  2.7217e+01 5.21695390
>> # of obs: 880, groups: Plot3, 20; Plot1:Plot2, 8; Plot1, 4
>>
>> Fixed effects:
>> Estimate Std. Error  DF t value Pr(>|t|)
>> (Intercept)  2.354550.24871 878  9.4671   <2e-16 ***
>> XvarFactor2  0.390910.35173 878  1.1114   0.2667
>> ---
>> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
>>
>> Look the wrong P value, I know that it is wrong because the DF  
>> used. But, In
>> this case, the result is not correct. Dont have any difference of  
>> the result
>> using random effects with lmer and using a simple analyses with lm.
>
> You are assuming that there is a correct value of the denominator
> degrees of freedom.  I don't believe there is.  The statistic that is
> quoted there doesn't have exactly an F distribution so there is no
> correct degrees of freedom.
>
> One thing you can do with lmer is to form a Markov Chain Monte Carlo
> sample from the posterior distribution of the parameters so you can
> check to see whether the value of zero is in the middle of the
> distribution of XvarFactor2 or not.
>
> It would be possible for me to recreate in lmer the rules used in lme
> for calculating denominator degrees of freedom associated with terms
> of the random effects.  However, the class of models fit by lmer is
> larger than the class of models fit by lme (at least as far as the
> structure of the random-effects terms goe

[R] ESS and Emacs

2005-12-29 Thread Mark Leeds
I have been using the document written by John Fox titled
Sn Introduction to ESS + XEmacs for Windows
Users of R.
 
It's a very nice document and
I went through it carefully but I got
an error when I finished it and launched XEmacs.
 
The error is "cannot open load file : ess-site".
 
So, I did more investigation
and it seems like there is a folder
 
Program Files/XEmacs/xemacs-packages/etc/ess
 
that should have been created when
I downloaded XEmacs but it doesn't exist. 
There is a folder called efs ( rather than ess ) but
I don't think that's it.
 
I tried installing XEmacs again but the
same thing happened.
 
Does anyone know anything about this ?
I am using Windows NT.
Thanks.
 
 Mark


**
This email and any files transmitted with it are confidentia...{{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


Re: [R] importing shapefiles into spatstat

2005-12-29 Thread Mulholland, Tom
You might also consider looking at the R-sig-Geo list which has lots of 
discussion about issues relating to file formats and the best ways to get data 
in and out the various packages that are used.

Tom

> -Original Message-
> From: [EMAIL PROTECTED]
> [mailto:[EMAIL PROTECTED] Behalf Of Charlotte Reemts
> Sent: Friday, 30 December 2005 3:59 AM
> To: r-help@stat.math.ethz.ch
> Subject: [R] importing shapefiles into spatstat
> 
> 
> Dear R users,
> I am using spatstat to analyze point patterns (tree 
> locations).  I would
> like to import the shapefile with the study area polygons 
> (six total) into R
> and use it to create the window for the spatstat analysis.  I 
> do not simply
> want to use a rectangle because the study areas spread out 
> over 4 ha.
> Any suggestions would be greatly appreciated.
> Thanks,
> Charlotte Reemts
> 
> 
> Charlotte Reemts
> Vegetation Ecologist
> The Nature Conservancy--Fort Hood (TX) Project
> P.O. Box 5190
> Fort Hood, TX 76544-0190
> 254-286-6745
> fax: 254-288-5039
> [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
>

__
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


Re: [R] use of tapply?

2005-12-29 Thread François Pinard
[tom wright]

> I'm still learning how to program with R and I was hoping someone 
> could take the time to show me how I can rewrite this code?

I'll try! :-)

>data.intersects<-data.frame(
>x=c(0.230,0.411,0.477,0.241,0.552,0.230),
>y=c(0.119,0.515,0.261,0.431,0.304,0.389),
>angle=vector(length=6),
>length=vector(length=6),
>row.names=c('tbr','trg','dbr','dbg','pbr','pbg'))


>calcDist<-function(x,y){
>#calcualates distance from origin (C)
>origin<-data.frame(x=0.34,y=0.36)
>dx<-origin$x-x
>dy<-origin$y-y

>length<-sqrt(dx^2+dy^2)
>angle<-asin(dy/length)
>return(list('length'=length,'angle'=angle))
>}

>for(iLoc in 1:length(data.intersects[,1])){
>result<-calcDist(data.intersects[iLoc,]$x,data.intersects[iLoc,]$y)
>data.intersects[iLoc,]$angle<-result$angle
>data.intersects[iLoc,]$length<-result$length
>}

Using `di' instead of `data.intersects' for short:

di <- data.frame(x=c(0.230, 0.411, 0.477, 0.241, 0.552, 0.230),
 y=c(0.119, 0.515, 0.261, 0.431, 0.304, 0.389),
 row.names=c('tbr', 'trg', 'dbr', 'dbg', 'pbr', 'pbg'))
di.c <- with(di, data.frame(x=x-0.34,  y=y-0.36))
di$length <- with(di.c, sqrt(x^2 + y^2))
di$angle <- with(di.c, atan2(y, x))

-- 
François Pinard   http://pinard.progiciels-bpi.ca

__
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


Re: [R] How to fit all points into plot?

2005-12-29 Thread Duncan Murdoch
On 12/29/2005 5:33 PM, Sean Davis wrote:
> 
> 
> On 12/29/05 5:19 PM, "Martin Lam" <[EMAIL PROTECTED]> wrote:
> 
> 
>>Hi,
>>
>>I have a problem when I want to add new points (or a
>>new line) to the graph. Some points (or parts of the
>>line) are not shown on the graph because they lie
>>beyond the scale of the axis. Is there a way to
>>overcome this so all points (or the entire line) are
>>shown on the graph? Here's an example of my problem:
>>
>>colors = c("red", "blue")
>>
>>plot(x=rnorm(100,0,1), y=runif(100), type="p", pch=21,
>>col = colors[1])
>>
>># if I add these points not all of them are shown on
>>the graph
>>lines(x=rnorm(100,3,1), y=runif(100), type="p",
>>pch=24, col = colors[2])
> 
> 
> You will need to use something like max and min for ALL of your x and y
> values and then set your own ylim and xlim for the plot.

This is the right idea, but a quicker way is to use the range() function:

> 
>  xold <- rnorm(100,0,1)
>  yold <- runif(100)
>  xnew <- rnorm(100,3,1)
>  ynew <- runif(100)
> 
>  plot(xold,yold,xlim=c(min(c(xold,xnew)),max(c(xold,xnew))),
   xlim=range(c(xold,xnew)),

>   ylim=c(min(c(yold,ynew)),max(c(yold,ynew
 ylim=range(c(yold,ynew))
> 
>  lines(xold,yold,type='p')
> 
> This isn't tested, but I hope you get the idea.
> 
> Sean
> 
> __
> 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

__
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


Re: [R] Which cluster function can be used to cluster a correlaiton matrix?

2005-12-29 Thread Earl F. Glynn
"Martin Maechler" <[EMAIL PROTECTED]> wrote in message
news:[EMAIL PROTECTED]

> The clue is to transform correlations to dissimilarities.
. . .
> and then use  hclust(), agnes(), pam(), [the latter two from
> package 'cluster'], ...
> with 'Dx' as dissimilarity

Perhaps this TechNote may be of interest:

Correlation "Distances" and Hierarchical Clustering
http://research.stowers-institute.org/efg/R/Visualization/cor-cluster/index.htm

efg

__
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


Re: [R] How to fit all points into plot?

2005-12-29 Thread Sean Davis



On 12/29/05 5:19 PM, "Martin Lam" <[EMAIL PROTECTED]> wrote:

> Hi,
> 
> I have a problem when I want to add new points (or a
> new line) to the graph. Some points (or parts of the
> line) are not shown on the graph because they lie
> beyond the scale of the axis. Is there a way to
> overcome this so all points (or the entire line) are
> shown on the graph? Here's an example of my problem:
> 
> colors = c("red", "blue")
> 
> plot(x=rnorm(100,0,1), y=runif(100), type="p", pch=21,
> col = colors[1])
> 
> # if I add these points not all of them are shown on
> the graph
> lines(x=rnorm(100,3,1), y=runif(100), type="p",
> pch=24, col = colors[2])

You will need to use something like max and min for ALL of your x and y
values and then set your own ylim and xlim for the plot.

 xold <- rnorm(100,0,1)
 yold <- runif(100)
 xnew <- rnorm(100,3,1)
 ynew <- runif(100)

 plot(xold,yold,xlim=c(min(c(xold,xnew)),max(c(xold,xnew))),
  ylim=c(min(c(yold,ynew)),max(c(yold,ynew

 lines(xold,yold,type='p')

This isn't tested, but I hope you get the idea.

Sean

__
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


[R] How to fit all points into plot?

2005-12-29 Thread Martin Lam
Hi,

I have a problem when I want to add new points (or a
new line) to the graph. Some points (or parts of the
line) are not shown on the graph because they lie
beyond the scale of the axis. Is there a way to
overcome this so all points (or the entire line) are
shown on the graph? Here's an example of my problem:

colors = c("red", "blue")

plot(x=rnorm(100,0,1), y=runif(100), type="p", pch=21,
col = colors[1])

# if I add these points not all of them are shown on
the graph
lines(x=rnorm(100,3,1), y=runif(100), type="p",
pch=24, col = colors[2])

Thanks in advance,

Martin

__
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


Re: [R] Axis/Ticks/Scale

2005-12-29 Thread Marc Schwartz (via MN)
On Thu, 2005-12-29 at 22:06 +0100, Martin Maechler wrote:
> > "Marc" == Marc Schwartz (via MN) <[EMAIL PROTECTED]>
> > on Wed, 28 Dec 2005 15:46:37 -0600 writes:
> 
> Marc> On Wed, 2005-12-28 at 20:15 +,
> Marc> [EMAIL PROTECTED] wrote:
> >> Dear All,
> >> 
> >> Apologies for this simple question and thanks in advance
> >> for any help given.
> >> 
> >> Suppose I wanted to plot 1 million observations and
> >> produce the command
> >> 
> >> plot(rnorm(100))
> >> 
> >> The labels of the xaxis are 0, e+00 2 e+05 etc. These are
> >> clearly not very attractive (The plots are for a
> >> PhD. thesis).
> >> 
> >> I would like the axes to be 0,2,4,6,8,10 with a *10^5 on
> >> the right hand side.
> >> 
> >> Is there a simple command for this?
> >> 
> >> Best Wishes
> >> 
> >> Roger
> 
> 
> Marc> See ?plotmath for some additional examples and there
> Marc> are some others in the list archives.
> 
> Yes, I think this one is there too:
> It has the "* 10^k" after each number;
> the nice thing about it is that it works for all kind of data
> -- and of course, in principle it could be built into R ...
> 
> 
> 
> ###- Do "a 10^k" labeling instead of "a e" ---
> 
> axTexpr <- function(side, at = axTicks(side, axp=axp, usr=usr, log=log),
> axp = NULL, usr = NULL, log = NULL)
> {
> ## Purpose: Do "a 10^k" labeling instead of "a e"
> ##  this auxiliary should return 'at' and 'label' (expression)
> ## --
> ## Arguments: as for axTicks()
> ## --
> ## Author: Martin Maechler, Date:  7 May 2004, 18:01
> eT <- floor(log10(abs(at)))# at == 0 case is dealt with below
> mT <- at / 10^eT
> ss <- lapply(seq(along = at),
>  function(i) if(at[i] == 0) quote(0) else
>  substitute(A %*% 10^E, list(A=mT[i], E=eT[i])))
> do.call("expression", ss)
> }
> 
> 
> x <- 1e7*(-10:50)
> y <- dnorm(x, m=10e7, s=20e7)
> plot(x,y)
> ## ^^ not so nice; ok, try 
> 
> par(mar=.1+c(5,5,4,1))## << For the horizontal y-axis labels, need more space
> plot(x,y, axes= FALSE, frame=TRUE)
> aX <- axTicks(1); axis(1, at=aX, label= axTexpr(1, aX))
> if(FALSE) # rather the next one
> { aY <- axTicks(2); axis(2, at=aY, label= axTexpr(2, aY))}
> ## or rather (horizontal labels on y-axis):
> aY <- axTicks(2); axis(2, at=aY, label= axTexpr(2, aY), las=2)


Nice Martin!

I do like that.  I also like the handling of zero, which I realized
after sending my reply, thus should have used:

  x <- rnorm(100)
  plot(x, xaxt = "n")
  x.at <- seq(0, 10, 2) * 10 ^ 5

  # Handle the zero here this time
  x.lab <- parse(text = paste(seq(2, 10, 2), "%*% 10^5"))
  axis(1, at = x.at, labels = c(0, x.lab))


BTW, on your approach, it was here:

  http://finzi.psych.upenn.edu/R/Rhelp02a/archive/36462.html

and more recently, here:

  http://finzi.psych.upenn.edu/R/Rhelp02a/archive/57011.html

:-)

Best regards and Happy New Year,

Marc

__
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


Re: [R] importing shapefiles into spatstat

2005-12-29 Thread Don MacQueen
And I see that I am probably already be out of date! There is a nice 
article in R News Volume 5/2, November 2005, available from CRAN.
-Don

At 1:41 PM -0800 12/29/05, Don MacQueen wrote:
>Go to CRAN, to the "Packages" page, do a simple search on the text "shape".
>It will quickly lead you to
>maptools   tools for reading and handling shapefiles
>  From there, I guess you'll have to extract the polygons from the
>structure that is returned by the function(s) in the maptools package.
>-Don
>
>
>
>At 1:59 PM -0600 12/29/05, Charlotte Reemts wrote:
>>Dear R users,
>>I am using spatstat to analyze point patterns (tree locations).  I would
>>like to import the shapefile with the study area polygons (six total) into R
>>and use it to create the window for the spatstat analysis.  I do not simply
>>want to use a rectangle because the study areas spread out over 4 ha.
>>Any suggestions would be greatly appreciated.
>>Thanks,
>>Charlotte Reemts
>>
>>
>>Charlotte Reemts
>>Vegetation Ecologist
>>The Nature Conservancy--Fort Hood (TX) Project
>>P.O. Box 5190
>>Fort Hood, TX 76544-0190
>>254-286-6745
>>fax: 254-288-5039
>>[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
>
>
>--
>--
>Don MacQueen
>Environmental Protection Department
>Lawrence Livermore National Laboratory
>Livermore, CA, USA
>
>__
>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


-- 
--
Don MacQueen
Environmental Protection Department
Lawrence Livermore National Laboratory
Livermore, CA, USA

__
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


Re: [R] importing shapefiles into spatstat

2005-12-29 Thread Don MacQueen
Go to CRAN, to the "Packages" page, do a simple search on the text "shape".
It will quickly lead you to
   maptools tools for reading and handling shapefiles
 From there, I guess you'll have to extract the polygons from the 
structure that is returned by the function(s) in the maptools package.
-Don



At 1:59 PM -0600 12/29/05, Charlotte Reemts wrote:
>Dear R users,
>I am using spatstat to analyze point patterns (tree locations).  I would
>like to import the shapefile with the study area polygons (six total) into R
>and use it to create the window for the spatstat analysis.  I do not simply
>want to use a rectangle because the study areas spread out over 4 ha.
>Any suggestions would be greatly appreciated.
>Thanks,
>Charlotte Reemts
>
>
>Charlotte Reemts
>Vegetation Ecologist
>The Nature Conservancy--Fort Hood (TX) Project
>P.O. Box 5190
>Fort Hood, TX 76544-0190
>254-286-6745
>fax: 254-288-5039
>[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


-- 
--
Don MacQueen
Environmental Protection Department
Lawrence Livermore National Laboratory
Livermore, CA, USA

__
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


Re: [R] trouble with S4 methods for group "Summary"

2005-12-29 Thread Martin Maechler
Yes, setting 'Summary'  S4 group methods is a bit painful,
because the S3 generic starts with "...".

In the 'Matrix' CRAN package,
we do the following  {thanks to hints by John Chambers IIRC}:

Our AllGeneric.R file
(https://svn.R-project.org/R-packages/trunk/Matrix/R/AllGeneric.R)
ends with

### Group Generics 
## The following are **WORKAROUND** s currently needed for all non-Primitives:

##  "Math"
setGeneric("log", group="Math")
setGeneric("gamma", group="Math")
setGeneric("lgamma", group="Math")

## "Math2"
setGeneric("round",  group="Math2")
setGeneric("signif", group="Math2")

## "Summary" --- this needs some hoop jumping that may become unnecessary
##   in a future version of R (>= 2.3.x):

.max_def <- function(x, ..., na.rm = FALSE) base::max(x, ..., na.rm = na.rm)
.min_def <- function(x, ..., na.rm = FALSE) base::min(x, ..., na.rm = na.rm)
.range_def <- function(x, ..., na.rm = FALSE) base::range(x, ..., na.rm = na.rm)
.prod_def <- function(x, ..., na.rm = FALSE) base::prod(x, ..., na.rm = na.rm)
.sum_def <- function(x, ..., na.rm = FALSE) base::sum(x, ..., na.rm = na.rm)
.any_def <- function(x, ..., na.rm = FALSE) base::any(x, ..., na.rm = na.rm)
.all_def <- function(x, ..., na.rm = FALSE) base::all(x, ..., na.rm = na.rm)

setGeneric("max", function(x, ..., na.rm = FALSE) standardGeneric("max"),
   useAsDefault = .max_def, group = "Summary")
setGeneric("min", function(x, ..., na.rm = FALSE) standardGeneric("min"),
   useAsDefault = .min_def, group="Summary")
setGeneric("range", function(x, ..., na.rm = FALSE) standardGeneric("range"),
   useAsDefault = .range_def, group="Summary")
setGeneric("prod", function(x, ..., na.rm = FALSE) standardGeneric("prod"),
   useAsDefault = .prod_def, group="Summary")
setGeneric("sum", function(x, ..., na.rm = FALSE) standardGeneric("sum"),
   useAsDefault = .sum_def, group="Summary")
setGeneric("any", function(x, ..., na.rm = FALSE) standardGeneric("any"),
   useAsDefault = .any_def, group="Summary")
setGeneric("all", function(x, ..., na.rm = FALSE) standardGeneric("all"),
   useAsDefault = .all_def, group="Summary")

##-

and then in dMatrix.R we have

## This needs extra work in ./AllGeneric.R :
setMethod("Summary", signature(x = "dMatrix", na.rm = "ANY"),
  function(x, ..., na.rm) callGeneric([EMAIL PROTECTED], ..., na.rm = 
na.rm))


I think you can safely follow this recipe;

Regards,
Martin Maechler, ETH Zurich


> "Parlamis" == Parlamis Franklin <[EMAIL PROTECTED]>
> on Wed, 28 Dec 2005 19:52:00 -1000 writes:

Parlamis> Hello.  This question concerns the Methods
Parlamis> package.  I have created a new class and am trying
Parlamis> to set a method for it for S4 group generic
Parlamis> "Summary".  I have run into some signature
Parlamis> problems.  An example:

>> setClass("track", representation(x="numeric",
>> y="character"))
Parlamis>   [1] "track"
>> setGeneric("max", group="Summary")
Parlamis>   [1] "max"
>> setMethod("Summary", signature(x="track"), function(x,
>> ..., na.rm)
Parlamis> callGeneric([EMAIL PROTECTED], ..., na.rm)) [1] "Summary"
>> dd<-new("track", x=c(1,2), y="abc") max(dd)
Parlamis>   [1] -Inf Warning message: no finite arguments to
Parlamis> max; returning -Inf
>> showMethods("max")

Parlamis>   Function "max": na.rm = "ANY" na.rm = "track"
Parlamis> na.rm = "missing" (inherited from na.rm = "ANY")

Parlamis> As you can see from the above, the method I tried
Parlamis> to set for "max" (via "Summary") was defined for
Parlamis> the formal argument "na.rm" not "x".  This makes
Parlamis> sense because the standardGeneric created for max
Parlamis> only allows methods to be defined for argument
Parlamis> "na.rm"

>> max
Parlamis>   standardGeneric for "max" defined from package
Parlamis> "base" belonging to group(s): Summary

Parlamis>   function (..., na.rm = FALSE)
Parlamis> standardGeneric("max") 
Parlamis> Methods may be defined for arguments: na.rm

Parlamis> However, group "Summary" purports to allow you to
Parlamis> define methods for arguments "x" and "na.rm".

>> Summary
Parlamis>   groupGenericFunction for "Summary" defined from
Parlamis> package "base"

Parlamis>   function (x, ..., na.rm = FALSE) stop("function
Parlamis> 'Summary' is a group generic; do not call it
Parlamis> directly", domain = NA) 
Parlamis> Methods may be defined for arguments: x, na.rm

Parlamis> How does this work?  Can someone point me to where
Parlamis> I am going wrong, and explain how to define S4
Parlamis> methods for group "Summary" for argument "x"?
Parlamis> Perhaps I need to do more in my "setGeneric" call?
Parlamis> Thanks in advance.

Parlamis> __
Parlamis> R-help@stat.math.et

Re: [R] Axis/Ticks/Scale

2005-12-29 Thread Martin Maechler
> "Marc" == Marc Schwartz (via MN) <[EMAIL PROTECTED]>
> on Wed, 28 Dec 2005 15:46:37 -0600 writes:

Marc> On Wed, 2005-12-28 at 20:15 +,
Marc> [EMAIL PROTECTED] wrote:
>> Dear All,
>> 
>> Apologies for this simple question and thanks in advance
>> for any help given.
>> 
>> Suppose I wanted to plot 1 million observations and
>> produce the command
>> 
>> plot(rnorm(100))
>> 
>> The labels of the xaxis are 0, e+00 2 e+05 etc. These are
>> clearly not very attractive (The plots are for a
>> PhD. thesis).
>> 
>> I would like the axes to be 0,2,4,6,8,10 with a *10^5 on
>> the right hand side.
>> 
>> Is there a simple command for this?
>> 
>> Best Wishes
>> 
>> Roger


Marc> See ?plotmath for some additional examples and there
Marc> are some others in the list archives.

Yes, I think this one is there too:
It has the "* 10^k" after each number;
the nice thing about it is that it works for all kind of data
-- and of course, in principle it could be built into R ...



###- Do "a 10^k" labeling instead of "a e" ---

axTexpr <- function(side, at = axTicks(side, axp=axp, usr=usr, log=log),
axp = NULL, usr = NULL, log = NULL)
{
## Purpose: Do "a 10^k" labeling instead of "a e"
##this auxiliary should return 'at' and 'label' (expression)
## --
## Arguments: as for axTicks()
## --
## Author: Martin Maechler, Date:  7 May 2004, 18:01
eT <- floor(log10(abs(at)))# at == 0 case is dealt with below
mT <- at / 10^eT
ss <- lapply(seq(along = at),
 function(i) if(at[i] == 0) quote(0) else
 substitute(A %*% 10^E, list(A=mT[i], E=eT[i])))
do.call("expression", ss)
}


x <- 1e7*(-10:50)
y <- dnorm(x, m=10e7, s=20e7)
plot(x,y)
## ^^ not so nice; ok, try 

par(mar=.1+c(5,5,4,1))## << For the horizontal y-axis labels, need more space
plot(x,y, axes= FALSE, frame=TRUE)
aX <- axTicks(1); axis(1, at=aX, label= axTexpr(1, aX))
if(FALSE) # rather the next one
{ aY <- axTicks(2); axis(2, at=aY, label= axTexpr(2, aY))}
## or rather (horizontal labels on y-axis):
aY <- axTicks(2); axis(2, at=aY, label= axTexpr(2, aY), las=2)

__
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


Re: [R] Problem Reading SPlus Dump Into R - Spaces Embedded in Data

2005-12-29 Thread Peter Dalgaard
allan miller <[EMAIL PROTECTED]> writes:

> Hello,
> 
> I'm trying to source() an SPlus 6.x file created using dump(..., 
> oldStyle=T) into R (version 2.01) as using the following instructions:
> 
> > *If you have access to S-PLUS, it is usually more reliable to |dump| 
> > the object(s) in S-PLUS and |source| the dumpfile in R. For S-PLUS 5.x 
> > and 6.x you may need to use |dump(..., oldStyle=T)|, and to read in 
> > very large objects it may be preferable to use the dumpfile as a batch 
> > script rather than use the |source| function.*
> 
> (from "R Data Import/Export," pg. 15)
> 
> An example:
> 
>  > source("testdump")
> Error in parse(file, n, text, prompt) : syntax error on line 1895
> 
> where the data on line 1895 - and other lines causing this - have 
> embedded spaces, such as the following:
> 
> 
> [line 1895]  Johnson Partners LLC
> 
> 
> I can't seem to find any options for either the SPlus dump, or R 
> source(), that relate to this problem.  Any suggestions for how to 
> either dump or source files containing data with embedded spaces?

A bit more context might be helpful. What's in lines surrounding 1895?
Can you show a simple S-PLUS object displaying the behaviour? What
happens if you dput() the object? Will S-PLUS itself restore the file? 

-- 
   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


[R] use of tapply?

2005-12-29 Thread tom wright
I'm still learning how to program with R and I was hoping someone could
take the time to show me how I can rewrite this code?
Many thanks
Tom

data.intersects<-data.frame(
x=c(0.230,0.411,0.477,0.241,0.552,0.230),
y=c(0.119,0.515,0.261,0.431,0.304,0.389),
angle=vector(length=6),
length=vector(length=6),
row.names=c('tbr','trg','dbr','dbg','pbr','pbg'))


calcDist<-function(x,y){
#calcualates distance from origin (C)
origin<-data.frame(x=0.34,y=0.36)
dx<-origin$x-x
dy<-origin$y-y

length<-sqrt(dx^2+dy^2)
angle<-asin(dy/length)
return(list('length'=length,'angle'=angle))
}

for(iLoc in 1:length(data.intersects[,1])){
result<-calcDist(data.intersects[iLoc,]$x,data.intersects[iLoc,]$y)
data.intersects[iLoc,]$angle<-result$angle
data.intersects[iLoc,]$length<-result$length
}

__
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


Re: [R] How to plot curves with more than 8 colors

2005-12-29 Thread Earl F. Glynn
"Don MacQueen" <[EMAIL PROTECTED]> wrote in message
news:[EMAIL PROTECTED]
>
> I have found this little function useful when trying to choose colors:

FWIW:  To choose colors, this page
http://research.stowers-institute.org/efg/R/Color/Chart/index.htm
and PDF
http://research.stowers-institute.org/efg/R/Color/Chart/ColorChart.pdf
can be used to view the "named" colors in R, along with RGB information in
decimal or hex.

Not everyone can "see" the same colors.  In particular, about 7% of males
have some color blindness, especially red-green.  See some of the links
under "Color Blindness":
http://www.efg2.com/Lab/Library/Color/index.html

Not all devices have the same color "gamut".  That is, not all devices can
display the same ranges of colors:
http://www.efg2.com/Lab/Library/Color/Science.htm#gamuts
Matching color between the screen and printers can be a problem.

efg

__
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


Re: [R] use of predict() with confidence/prediction bands

2005-12-29 Thread Peter Dalgaard
Alan Arnholt <[EMAIL PROTECTED]> writes:

> To my understanding, a confidence interval typically covers a single
> valued parameter.  In contrast, a confidence band covers an entire line
> with a band.  In regression, it is quite common to construct confidence
> and prediction bands.  I have found that many people are connecting
> individual confidence/prediction interval values produced with
> predict(object,sd.fit=T,type="conf/pred") and calling the result a
> confidence/prediction band.  Since there is no specific probability
> statement that can be attached to these connected confidence/prediction
> intervals, this does not seem reasonable to me.  This is done, for
> example, in ISWR pg. 105, UsingR for Introductory Statistics pg 296, and
> Linear Models with R pg. 39 (Although in this instance the intervals are
> called 95% "pointwise" confidence bands versus simply 95% confidence
> bands.)  To make a confidence/prediction band, one should  construct
> simultaneous confidence/prediction intervals with say a Scheffe approach
> as mentioned in the S-PLUS Guide to statistics pg 274.  If these connected
> intervals were called pointwise confidence/prediction intervals with the
> understanding that have no particular probability interpretation, then
> they are useful in understanding where the line should fall.  However,
> they are not confidence/prediction bands as such, and I think it is
> misleading to name them so.  Should the intervals the authors in the
> three mentioned references construct not be called something similar
> to connected 95% pointwise confidence/prediction intervals versus 95%
> confidence/prediction bands?  Or, have I missed the boat?  Fire away...

You do have a point, of course. My take is that (a) they are bands and
(b) they have the property that for _each_ x they contain y(x) with
95% probability. So "95% pointwise confidence bands" is reasonably
warranted to my mind. ISwR could probably be more careful in making
the "pointwise" distinction, but I'd be afraid of confusing readers
who might well be at the level where the prime difficulty is grasping
the difference between prediction intervals and confidence intervals.

Global coverage, i.e., bands that contain the true line with 95%
probability, is quite a bit harder to obtain, especially in the
nonparametric regression extensions. Such bands end up being rather
wide, and some (I'm afraid I forgot who) have suggested just to use
the pointwise bands with the understanding that they cover, on
average, 95% of the true line.

-- 
   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


[R] importing shapefiles into spatstat

2005-12-29 Thread Charlotte Reemts
Dear R users,
I am using spatstat to analyze point patterns (tree locations).  I would
like to import the shapefile with the study area polygons (six total) into R
and use it to create the window for the spatstat analysis.  I do not simply
want to use a rectangle because the study areas spread out over 4 ha.
Any suggestions would be greatly appreciated.
Thanks,
Charlotte Reemts


Charlotte Reemts
Vegetation Ecologist
The Nature Conservancy--Fort Hood (TX) Project
P.O. Box 5190
Fort Hood, TX 76544-0190
254-286-6745
fax: 254-288-5039
[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


Re: [R] loop

2005-12-29 Thread Peter Dalgaard
Duncan Murdoch <[EMAIL PROTECTED]> writes:

> On 12/29/2005 7:41 AM, gynmeerut wrote:
> > Dear All,
> > 
> > I have to use loop over an array  so I am using following procedure
> >  
> > count<-1
> >  repeat{
> >  count<-count + 1
> >  c(g[count],1:i[count]) ->qw
> >  if(count>5)break
> >  }
> 
> We can't reproduce this, as we don't have g or i.  But the general 
> advice in a case like this is to simulate the loop by hand:  write down 
> what is in each of the variables, and walk through the loop.
> 
> I suspect your problem is in initializing count improperly, or in 
> putting the test in the wrong place.

I think not. The assignment to "qw" only depends on "count", so it's small
wonder that the end result is that of the last iteration. If the
intention was to append to "qw", then you need something like

  qw <- c(qw, g[count],1:i[count])
 

> Duncan Murdoch
> > 
> > as  a result qw is
> > [1]  0.9643836  1.000  2.000  3.000  4.000  5.000
> >  [7]  6.000  7.000  8.000  9.000 10.000 11.000
> > [13] 12.000 13.000 14.000 15.000 16.000 17.000
> > [19] 18.000 19.000 20.000 21.000 22.000 23.000
> > [25] 24.000 25.000 26.000 27.000 28.000 29.000
> > [31] 30.000 31.000 32.000 33.000 34.000 35.000
> > [37] 36.000 37.000 38.000 39.000
> > 
> > which is according to the last value of  count i.e 6. What is the problem 
> > in error
> >  
> > 
> > but I expect the result as
> > 
> >  [1]  0.8328767  0.2410959  0.5315068  0.1424658  0.2520548  0.9643836
> >  [7]  6.000  7.000  8.000  0.2410959  1.000  2.000
> > [13]  3.000  4.000  5.000  0.5315068  1.000  2.000
> > [19]  3.000  4.000  5.000  6.000  7.000  8.000
> > [25]  9.000 10.000 11.000 12.000 13.000  0.1424658
> > [31]  1.000  0.000  0.2520548  1.000  2.000  3.000
> > [37]  4.000  5.000  6.000  7.000  8.000  9.000
> > [43] 10.000 11.000 12.000 13.000 14.000 15.000
> > [49] 16.000 17.000 18.000 19.000 20.000 21.000
> > [55] 22.000 23.000  0.9643836  1.000  2.000  3.000
> > [61]  4.000  5.000  6.000  7.000  8.000  9.000
> > [67] 10.000 11.000 12.000 13.000 14.000 15.000
> > [73] 16.000 17.000 18.000 19.000 20.000 21.000
> > [79] 22.000 23.000 24.000 25.000 26.000 27.000
> > [85] 28.000 29.000 30.000 31.000 32.000 33.000
> > [91] 34.000 35.000 36.000 37.000 38.000 39.000
> > 
> > which is combination of 6 vectors of lengths {9,6,14,3,24,40).
> > 
> > 
> > Please tell me the error or alternate procedure for looping.
> > 
> > regards,
> > 
> > GS
> > 
> > __
> > 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
> 
> __
> 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
> 

-- 
   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


Re: [R] function cv.glm in library 'boot'

2005-12-29 Thread Prof Brian Ripley
Please do RTFM.  It uses the cost function given by its 'cost' argument.

Using your suggestion to choose the threshold is not honest (in the 
technical sense of the word).

On Thu, 29 Dec 2005, Jun Ding wrote:

> Hi, everyone,
> I have a question regarding function cv.glm in library
> 'boot'.
> Basically cv.glm can  calculate the estimated K-fold
> cross-validation prediction error for generalized
> linear models. My question is this: if I am fitting a
> logit model, what kind of threshold will it use to
> calculate the prediction error (saved in 'delta')? It
> will use 0.5 as the threshold or pick a threshold such
> that the prediction error reaches its minimum?
>
> Thank you!
>
> Best wishes,
> Jun
>
> __
> 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

PLEASE do take note of that, as you asked a question which would have 
been answered by doing the homework asked for there.

-- 
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


[R] Problem Reading SPlus Dump Into R - Spaces Embedded in Data

2005-12-29 Thread allan miller
Hello,

I'm trying to source() an SPlus 6.x file created using dump(..., 
oldStyle=T) into R (version 2.01) as using the following instructions:

> *If you have access to S-PLUS, it is usually more reliable to |dump| 
> the object(s) in S-PLUS and |source| the dumpfile in R. For S-PLUS 5.x 
> and 6.x you may need to use |dump(..., oldStyle=T)|, and to read in 
> very large objects it may be preferable to use the dumpfile as a batch 
> script rather than use the |source| function.*

(from "R Data Import/Export," pg. 15)

An example:

 > source("testdump")
Error in parse(file, n, text, prompt) : syntax error on line 1895

where the data on line 1895 - and other lines causing this - have 
embedded spaces, such as the following:


[line 1895]  Johnson Partners LLC


I can't seem to find any options for either the SPlus dump, or R 
source(), that relate to this problem.  Any suggestions for how to 
either dump or source files containing data with embedded spaces?

Thank You.


[[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


Re: [R] Show graph integrated to GUI

2005-12-29 Thread Peter Dalgaard
"Michael H. Prager" <[EMAIL PROTECTED]> writes:

> I'm seeking suggestions for a book on Tcl/Tk, as it can be used with R. 
> The book I bought for the purpose ("Effective Tcl/Tk Programming") seems 
> quite unsuitable. For example, it has no description of anything like 
> the slider control provided by tkscale().
> 
> I have located the Tcl and Tk documentation provided with the R 
> distribution. That has plenty of detail. Now I would like something with 
> more of an overview, including numerous examples and illustrations.

The one book on the subject that I have and like (out of a total of
two...) is Welch's "Practical Programming in Tcl and Tk". I have the
3rd ed., the 4th ed is said to be even better (ISBN: 0130385603).
There are others though; you could do worse than consulting
wiki.tcl.tk/57 and the newsgroup comp.lang.tcl and archives of same.

 
> Thanks,
> MHP
> 
> 
> on 12/29/2005 3:06 AM Prof Brian Ripley said the following:
> 
> >This can be done with tkrplot package, whose sole example is of a slider 
> >and a graph in the same window.
> >
> >(AFAICS the slider function in the TeachingDemos package shows no more 
> >concepts than the tkdensity demo which ships with R.)
> >
> >On Thu, 29 Dec 2005, Chihiro Kuraya wrote:
> >
> >  
> >
> >>Hi,
> >>
> >>Thank you for information..
> >>
> >>I tried the TechingDemos package.
> >>But it seems that the window which contains slider control
> >>and graph window is separated.
> >>
> >>What I want to do is show slider control and graph
> >>in one window simultaneiously.
> >>It is possible ?
> >>
> >>Chihiro Kuraya
> >>
> >>
> >>"Gregory Snow" <[EMAIL PROTECTED]> wrote:
> >>
> >>
> >>
> >>>The slider function in the TeachingDemos and relax packages (same
> >>>function is in both packages you can use either) provides a way to do
> >>>this using a Tk window.  There are also several functions in the
> >>>TechingDemos package that use a lower level interface to a Tk window
> >>>that you can look at their source code as an example to build your own
> >>>(examples include: vis.gamma, rotate.persp, and run.power.examp).
> >>>
> >>>Hope this helps,
> >>>  
> >>>
> >
> >  
> >
> -Original Message-
> From: [EMAIL PROTECTED]
> [mailto:[EMAIL PROTECTED] On Behalf Of Chihiro Kuraya
> Sent: Saturday, December 24, 2005 9:53 PM
> To: r-help@stat.math.ethz.ch
> Subject: [R] Show graph integrated to GUI
> 
> Hi all,
> 
> It it posssible to show graph which is integrated to some GUI
> (e.g. TclTk or R-wxPython).
> 
> I want to make an application by R,
> for example, like the following picture:
> http://www.natch.co.uk/downloads/SigJenny/SJnScreenShot.gif
> 
> 
> >
> >  
> >
> 
> -- 
> 
> Michael Prager, Ph.D.
> Population Dynamics Team, NMFS SE Fisheries Science Center
> NOAA Center for Coastal Fisheries and Habitat Research
> Beaufort, North Carolina  28516
> http://shrimp.ccfhrb.noaa.gov/~mprager/
> Opinions expressed are personal, not official.  No
> government endorsement of any product is made or implied.
> 
> __
> 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
> 

-- 
   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


Re: [R] Testing a linear hypothesis after maximum likelihood

2005-12-29 Thread Spencer Graves
  1.  I try to avoid dogmatism and use whatever seems sufficiently 
accurate for the intended purposes and easiest to explain to the 
intended audience.

  2.  I'm not aware of any package that will compute Wald tests from 
optim(...)$hessian, etc., so I write my own code when I want that.

  3.  Likelihood ratio tests are known to be more accurate than Wald 
tests.  Linear regression can be thought of as projection onto a 
subspace.  Nonlinear least squares and maximum likelihood more generally 
involve projection onto a nonlinear manifold.  It does this by creating 
local linear approximations.  There are two sources of error in this due 
to (1) intrinsic curvature of the manifold and (2) parameter effects 
curvature.  I mention this, because likelihood ratio procedures are 
distorted only by the intrinsic curvature, while Wald procedures are 
subject to both.  Moreover, in evaluating numerous published 
applications of nonlinear least squares, Bates and Watts found that the 
intrinsic curvature was never much worse than the parameter effects and 
was usually at least an order of magnitude smaller. See Bates and Watts 
(1988) Nonlinear Regression Analysis and Its Applications (Wiley) or 
Seber and Wild (1988) Nonlinear Regression (Wiley).

  Bottom line:  I routinely use Wald procedures to compute confidence 
intervals, because computing them by profiling log(likelihood ratio) is 
usually more work than I have time for.  However, for testing, when I 
have the time, I use likelihood ratio procedures.

  spencer graves

Peter Muhlberger wrote:
 > On 12/29/05 1:35 PM, "Spencer Graves" <[EMAIL PROTECTED]> wrote:
 >
 >
 >> I think the question was appropriate for this list.  If you want to
 >>do a Wald test, you might consider asking "optim" for "hessian=TRUE".
 >>If the function that "optim" minimizes is (-log(likelihood)), then the
 >>optional component "hessian" of the output of optim should be the
 >>observed information matrix.  An inverse of that should then estimate
 >>the parameter covariance matrix.  I often use that when "nls" dies on
 >>me, because "optim" will give me an answer.  If the hessian is singular,
 >>I can sometimes diagnose the problem by looking at eigenvalues and
 >>eigenvectors of the hessian.
 >
 >
 > Niffty, thanks again!  Do you construct your own wald tests out of 
matrixes
 > or use something packaged?  Or do you just avoid wald tests at all 
costs :)
 > ?
 >
 > Peter
 >
Spencer Graves wrote:

>   I think the question was appropriate for this list.  If you want 
> to do a Wald test, you might consider asking "optim" for "hessian=TRUE". 
> If the function that "optim" minimizes is (-log(likelihood)), then the 
> optional component "hessian" of the output of optim should be the 
> observed information matrix.  An inverse of that should then estimate 
> the parameter covariance matrix.  I often use that when "nls" dies on 
> me, because "optim" will give me an answer.  If the hessian is singular, 
> I can sometimes diagnose the problem by looking at eigenvalues and 
> eigenvectors of the hessian.
> 
>   hope this helps.
>   spencer graves
> 
> 
> On 12/29/05 7:04 AM, "Spencer Graves" <[EMAIL PROTECTED]> wrote:
> 
> 
>  >>  Why can't you use a likelihood ratio?  I would write two slightly
>  >> different functions, the second of which would use the linear 
> constraint
>  >> to eliminate one of the coefficients.  Then I'd refer 2*log(likelihood
>  >> ratio) to chi-square(1).  If I had some question about the chi-square
>  >> approximation to the distribution of that 2*log(likelihood ratio)
>  >> statistic, I'm use some kind of Monte Carlo, e.g., MCMC.
>  >>
> 
> 
> Neat solution, thanks!  I didn't see that, having focused my attention on
> finding some way to do a Wald test.  I think I was so focused because I
> thought it would be good to have some way of testing hypotheses w/o having
> to rerun my model every time.
> 
> 
>  >>  If you'd like more help from this listserve, PLEASE do read the
>  >> posting guide! "www.R-project.org/posting-guide.html".  Anecdotal
>  >> evidence suggests that posts that follow more closely the 
> suggestions in
>  >> that guide tend to get more useful replies quicker.
> 
> 
> Ok, I guess you're hinting that I'm violating the 'do your homework' norm.
> I'm not a statistician (I'm a social scientist) & was thinking about
> alternatives to the likelihood ratio test, so the self-evident solution you
> mention above didn't occur to me.  I did spend a long time trying to figure
> out whether there were facilities for Wald tests and whether they might 
> work
> w/ ML output.  It wasn't clear what would work & it would have taken even
> more time to try some alternatives out, so I thought I'd just ask the
> list--surely people have tests they typically run after ML.
> 
> In hindsight, I guess the question as asked was rather dumb, so my
> apologies.  Perhaps I should have asked if anyone

[R] use of predict() with confidence/prediction bands

2005-12-29 Thread Alan Arnholt


To my understanding, a confidence interval typically covers a single
valued parameter.  In contrast, a confidence band covers an entire line
with a band.  In regression, it is quite common to construct confidence
and prediction bands.  I have found that many people are connecting
individual confidence/prediction interval values produced with
predict(object,sd.fit=T,type="conf/pred") and calling the result a
confidence/prediction band.  Since there is no specific probability
statement that can be attached to these connected confidence/prediction
intervals, this does not seem reasonable to me.  This is done, for
example, in ISWR pg. 105, UsingR for Introductory Statistics pg 296, and
Linear Models with R pg. 39 (Although in this instance the intervals are
called 95% "pointwise" confidence bands versus simply 95% confidence
bands.)  To make a confidence/prediction band, one should  construct
simultaneous confidence/prediction intervals with say a Scheffe approach
as mentioned in the S-PLUS Guide to statistics pg 274.  If these connected
intervals were called pointwise confidence/prediction intervals with the
understanding that have no particular probability interpretation, then
they are useful in understanding where the line should fall.  However,
they are not confidence/prediction bands as such, and I think it is
misleading to name them so.  Should the intervals the authors in the
three mentioned references construct not be called something similar
to connected 95% pointwise confidence/prediction intervals versus 95%
confidence/prediction bands?  Or, have I missed the boat?  Fire away...

Alan T. Arnholt
Associate Professor
Dept. of Mathematical Sciences
Appalachian State University

__
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


[R] function cv.glm in library 'boot'

2005-12-29 Thread Jun Ding
Hi, everyone, 
I have a question regarding function cv.glm in library
'boot'. 
Basically cv.glm can  calculate the estimated K-fold
cross-validation prediction error for generalized
linear models. My question is this: if I am fitting a
logit model, what kind of threshold will it use to
calculate the prediction error (saved in 'delta')? It
will use 0.5 as the threshold or pick a threshold such
that the prediction error reaches its minimum?

Thank you!

Best wishes,
Jun

__
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


Re: [R] Show graph integrated to GUI

2005-12-29 Thread Michael H. Prager
I'm seeking suggestions for a book on Tcl/Tk, as it can be used with R. 
The book I bought for the purpose ("Effective Tcl/Tk Programming") seems 
quite unsuitable. For example, it has no description of anything like 
the slider control provided by tkscale().

I have located the Tcl and Tk documentation provided with the R 
distribution. That has plenty of detail. Now I would like something with 
more of an overview, including numerous examples and illustrations.

Thanks,
MHP


on 12/29/2005 3:06 AM Prof Brian Ripley said the following:

>This can be done with tkrplot package, whose sole example is of a slider 
>and a graph in the same window.
>
>(AFAICS the slider function in the TeachingDemos package shows no more 
>concepts than the tkdensity demo which ships with R.)
>
>On Thu, 29 Dec 2005, Chihiro Kuraya wrote:
>
>  
>
>>Hi,
>>
>>Thank you for information..
>>
>>I tried the TechingDemos package.
>>But it seems that the window which contains slider control
>>and graph window is separated.
>>
>>What I want to do is show slider control and graph
>>in one window simultaneiously.
>>It is possible ?
>>
>>Chihiro Kuraya
>>
>>
>>"Gregory Snow" <[EMAIL PROTECTED]> wrote:
>>
>>
>>
>>>The slider function in the TeachingDemos and relax packages (same
>>>function is in both packages you can use either) provides a way to do
>>>this using a Tk window.  There are also several functions in the
>>>TechingDemos package that use a lower level interface to a Tk window
>>>that you can look at their source code as an example to build your own
>>>(examples include: vis.gamma, rotate.persp, and run.power.examp).
>>>
>>>Hope this helps,
>>>  
>>>
>
>  
>
-Original Message-
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] On Behalf Of Chihiro Kuraya
Sent: Saturday, December 24, 2005 9:53 PM
To: r-help@stat.math.ethz.ch
Subject: [R] Show graph integrated to GUI

Hi all,

It it posssible to show graph which is integrated to some GUI
(e.g. TclTk or R-wxPython).

I want to make an application by R,
for example, like the following picture:
http://www.natch.co.uk/downloads/SigJenny/SJnScreenShot.gif


>
>  
>

-- 

Michael Prager, Ph.D.
Population Dynamics Team, NMFS SE Fisheries Science Center
NOAA Center for Coastal Fisheries and Habitat Research
Beaufort, North Carolina  28516
http://shrimp.ccfhrb.noaa.gov/~mprager/
Opinions expressed are personal, not official.  No
government endorsement of any product is made or implied.

__
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


Re: [R] Testing a linear hypothesis after maximum likelihood

2005-12-29 Thread Spencer Graves
  I think the question was appropriate for this list.  If you want to 
do a Wald test, you might consider asking "optim" for "hessian=TRUE". 
If the function that "optim" minimizes is (-log(likelihood)), then the 
optional component "hessian" of the output of optim should be the 
observed information matrix.  An inverse of that should then estimate 
the parameter covariance matrix.  I often use that when "nls" dies on 
me, because "optim" will give me an answer.  If the hessian is singular, 
I can sometimes diagnose the problem by looking at eigenvalues and 
eigenvectors of the hessian.

  hope this helps.
  spencer graves


On 12/29/05 7:04 AM, "Spencer Graves" <[EMAIL PROTECTED]> wrote:


 >>  Why can't you use a likelihood ratio?  I would write two slightly
 >> different functions, the second of which would use the linear constraint
 >> to eliminate one of the coefficients.  Then I'd refer 2*log(likelihood
 >> ratio) to chi-square(1).  If I had some question about the chi-square
 >> approximation to the distribution of that 2*log(likelihood ratio)
 >> statistic, I'm use some kind of Monte Carlo, e.g., MCMC.
 >>


Neat solution, thanks!  I didn't see that, having focused my attention on
finding some way to do a Wald test.  I think I was so focused because I
thought it would be good to have some way of testing hypotheses w/o having
to rerun my model every time.


 >>  If you'd like more help from this listserve, PLEASE do read the
 >> posting guide! "www.R-project.org/posting-guide.html".  Anecdotal
 >> evidence suggests that posts that follow more closely the suggestions in
 >> that guide tend to get more useful replies quicker.


Ok, I guess you're hinting that I'm violating the 'do your homework' norm.
I'm not a statistician (I'm a social scientist) & was thinking about
alternatives to the likelihood ratio test, so the self-evident solution you
mention above didn't occur to me.  I did spend a long time trying to figure
out whether there were facilities for Wald tests and whether they might work
w/ ML output.  It wasn't clear what would work & it would have taken even
more time to try some alternatives out, so I thought I'd just ask the
list--surely people have tests they typically run after ML.

In hindsight, I guess the question as asked was rather dumb, so my
apologies.  Perhaps I should have asked if anyone uses a built-in Wald
function after ML?  Or perhaps even that question is far too basic for a
list composed of such capable people.

Anyway, thanks for the insight!

Peter
#
  Why can't you use a likelihood ratio?  I would write two slightly
different functions, the second of which would use the linear constraint
to eliminate one of the coefficients.  Then I'd refer 2*log(likelihood
ratio) to chi-square(1).  If I had some question about the chi-square
approximation to the distribution of that 2*log(likelihood ratio)
statistic, I'm use some kind of Monte Carlo, e.g., MCMC.

  If you'd like more help from this listserve, PLEASE do read the
posting guide! "www.R-project.org/posting-guide.html".  Anecdotal
evidence suggests that posts that follow more closely the suggestions in
that guide tend to get more useful replies quicker.

  hope this helps.
  spencer graves


Peter Muhlberger wrote:

> I'd like to be able to test linear hypotheses after setting up and running a
> model using optim or perhaps nlm.  One hypothesis I need to test are that
> the average of several coefficients is less than zero, so I don't believe I
> can use the likelihood ratio test.
> 
> I can't seem to find a provision anywhere for testing linear combinations of
> coefficients after max. likelihood.
> 
> Cheers & happy holidays,
> 
> Peter
> 
> __
> 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

-- 
Spencer Graves, PhD
Senior Development Engineer
PDF Solutions, Inc.
333 West San Carlos Street Suite 700
San Jose, CA 95110, USA

[EMAIL PROTECTED]
www.pdf.com 
Tel:  408-938-4420
Fax: 408-280-7915

__
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


Re: [R] update to posting guide: use 'sessionInfo()' instead of 'version'

2005-12-29 Thread Seth Falcon
I think using sessionInfo() instead of version is a good idea.

On 29 Dec 2005, [EMAIL PROTECTED] wrote:
> [Note that sessionInfo() currently does not report all the
> information that 'version' does (it omits at least "Status" and "svn
> rev").  R-core members are aware of this -- whether or not they
> change this is up to them.]

Another thing not currently reported by sessionInfo() is a list of
loaded name spaces.  As more packages use name spaces and importing
dependencies instead of attaching them, this will become useful debug
info.

Loaded name spaces can be listed using loadedNamespaces().

+ seth

__
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


[R] addition of error terms

2005-12-29 Thread Manuel Gutierrez
Are there any functions to do error propagation in R?
 I have done a search with little success. Any
pointers to read about this topic are greatly
welcomed.
My specific problem is: I use a linear model (lm) to
predict the biomass of an individual in a population,
then I add up the biomass of all individuals to
calculate the total mass of the population. I want to
calculate the error in the final estimate of the total
biomass.
Also, I have another case where I use a nls model
instead of the lm, but in nls se.fit is currently
ignored as stated in the help page. Is there an
alternative way to calculate the errors in the
predictions from nls?
Thanks,
Manuel

__
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


[R] Glimmix and glm

2005-12-29 Thread Antonio_Paredes
Hello.

Some months age an e-mail was posted in which a comparison between Glimmix 
and glm was discussed. I have not been able to find that e-mail on the R 
archive. Does anyone recall the date of the above e-mail?

Thank you very much.

***
Antonio Paredes
USDA- Center for Veterinary Biologics
Biometrics Unit
510 South 17th Street, Suite 104
Ames, IA 50010
(515) 232-5785

 
[[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


[R] error propagation

2005-12-29 Thread Manuel Gutierrez
Are there any functions to do error propagation in R?
 I have done a search with little success. Any
pointers to read about this topic are greatly
welcomed.
My specific problem is: I use a linear model (lm) to
predict the biomass of an individual in a population,
then I add up the biomass of all individuals to
calculate the total mass of the population. I want to
calculate the error in the final estimate of the total
biomass.
Also, I have another case where I use a nls model
instead of the lm, but in nls se.fit is currently
ignored as stated in the help page. Is there an
alternative way to calculate the errors in the
predictions from nls?
Thanks,
Manuel

__
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


[R] update to posting guide: use 'sessionInfo()' instead of 'version'

2005-12-29 Thread Tony Plate
Some changes have been made to the posting guide, based on suggestions 
from various R-help contributors over the past year.

The most significant change is the recommendation to use 'sessionInfo()' 
  rather than 'version' when asking questions about unexpected behavior 
or bugs.  This change was made because 'sessionInfo()' reports the 
version and a list of packages currently attached.  As more and more 
packages become available, it becomes more likely that unexpected 
behavior is due to conflicts between packages, so this is relevant 
information.

[Note that sessionInfo() currently does not report all the information 
that 'version' does (it omits at least "Status" and "svn rev").  R-core 
members are aware of this -- whether or not they change this is up to them.]

-- Tony Plate

__
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


[R] Forward reference in Sweave

2005-12-29 Thread Dieter Menne
Dear Rweavers,

When generating reports with Sweave, I would like to quote some results in
the abstract (Something like "The treatment effect is 10 mbar, see page
100).

Currently, I used verbatimwrite and friends to write to multiple files to be
included,  but I wonder is there is a more elegant method for such
accumulated forward references in one file similar to toc creation.

Sorry if this is a latex-question, but it looks like most latexer are happy
with toc and index forward references. Maybe I used the wrong keywords for a
search.

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


Re: [R] Getting the same y-axis in a multivariate time series plot - plot(ts(...)); ylim does not do the trick

2005-12-29 Thread Gabor Grothendieck
If you use plot.zoo it will do it.  Just insert as.zoo(...) around
the ts object.

library(zoo)
plot(as.zoo(ts(cbind(rnorm(10), rnorm(10,mean=4,ylim=c(0,20))

On 12/29/05, Søren Højsgaard <[EMAIL PROTECTED]> wrote:
> I try to obtain the same y-axis for a 2-dim time series with
>
>  plot(ts(cbind(rnorm(10), rnorm(10,mean=4))),ylim=c(0,20))
>
> but that does not work. Looking in the code for plot.ts, the ylim-argument 
> seems to be taken care of, but not the way I expect. Can anyone help on this?
> Thanks
> Søren
>
> __
> 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
>

__
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


Re: [R] S4 classes: referencing slots with other slots

2005-12-29 Thread Gabor Grothendieck
Although this does not answer your question regarding how
to do it in S4, it is simple to do with the proto package.

This creates a fooWfcn proto object and then two child
proto objects of it.  Here foo inherits dat1 from fooWfcn
but foo2 has its own dat1 which overrides the dat1 in
its parent.

library(proto)
fooWfcn <- proto(dat1 = 1:10, dat2 = 10:20,
fn1 = function(., x) x - mean(.$dat1),
fn2 = function(., x) mean(.$dat2))

foo <- fooWfcn$proto()
foo$fn1(0) # -5.5

foo2 <- fooWfcn$proto(dat1 = 11:20)
foo2$fn1(0) # -15.5



On 12/29/05, A.J. Rossini <[EMAIL PROTECTED]> wrote:
> For those who suggest other ways to do this, I ALREADY HAVE ANOTHER
> DESIGN SOLUTION, DESCRIBED AT THE END.
>
> That being said, I want to know if it's possible to reference a slot
> in an S4 class from another slot, i.e. I'd like to have the "self.*"
> semantics of Python so that I can reuse a slot.  That is, for various
> reasons it would be nice to be able to do something like:
>
> setClass("fooWfcn",
> representation(dat1="vector",
>dat2="vector",
>fn1="function",
>fn2="function"),
> prototype=list(dat1=0:10,
>   dat2=10:20,
>   fn1=function(x) { return(x - mean(self.dat1)) },
>   fn2=function() { mean(self.dat2) }))
>
> and in the context of
>
> foo <- new("fooWfcn")
>
> have self.dat2 refer to [EMAIL PROTECTED], etc (I.e. in an instantiated 
> object,
> have the reference be within the object).
>
> One could easily imagine doing this on the fly, as well, during the "new" 
> call.
>
>  I've looked this up in a number of books, on-line talks/papers, etc,
> but havn't managed to find the right page describing it as possible or
> impossible.  I think it is the latter (impossible), since one could
> easily end up with a nasty self-referencing infinite loop (fn1
> refering to fn2 refering to fn1).  If it's the former, I'd be
> interested in knowing about it.
>
> Anyway, here's the DESIGN SOLUTION: write methods which do the
> appropriate "combination".
>
> Critique:
>
> pro - it works, it's simple, and I've already done it.
>
> con - for the problem I'm looking at, it's not quite so clean, adding
> one more layer of indirection that in Python or CLOS I'd not need,
> multiplied by a fair number of subclasses.
>
> best,
> -tony
>
> [EMAIL PROTECTED]
> Muttenz, Switzerland.
> "Commit early,commit often, and commit in a repository from which we can 
> easily
> roll-back your mistakes" (AJR, 4Jan05).
>
> __
> 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
>

__
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


Re: [R] How to plot curves with more than 8 colors

2005-12-29 Thread Don MacQueen

I have found this little function useful when trying to choose colors:

showcols <-  function (indx = 0:6)
{
 for (ii in unique(indx)) {
 is <- 100 * ii + 1:100
 if (min(is) > length(colors())) {
 cat("Maximum value of arg is", floor(length(colors())/100),
 "\n")
 return(NULL)
 }
 foo <- matrix(colors()[is], nrow = 10)
 par(mar = c(3, 3, 0.25, 0.25))
 plot(1:10, 1:10, type = "n", yaxt = "n", xlab = "", ylab = "")
 axis(2, at = 1:10, lab = 10:1)
 for (j in 1:10) {
 for (i in 1:10) {
 points(j, 11 - i, col = foo[i, j], pch = 16,
   cex = 4)
 text(j, 11 - i - 0.3, foo[i, j], cex = 0.8)
 }
 }
 if (length(indx) > 1 & ii < max(indx))
 readline(paste("Currently showing group", ii, "  CR to continue "))
 }
 invisible(foo)
}

Just type
showcols()
at the R prompt. Then pick any 20 you happen to think are distinguishable.

Based on my own experience, I would doubt that it is possible to find 
20 easily distinguishable colors.
But perhaps if you use both solid and dotted lines you can get 20 
distinguishable lines.

-Don

At 10:48 AM +0800 12/28/05, Vincent Deng wrote:
>Hi,
>
>Thanks for your kindly reply.
>I think maybe I didn't specify color codes properly. That is,the
>difference between each color is not sharp enough for me to identify
>them as different colors.
>
>So can you tell me about how to specify the color properly so that the
>difference among each color can be identified clearly?
>
>Thanks again and again ...
>
>On 12/27/05, Uwe Ligges <[EMAIL PROTECTED]> wrote:
>>  Vincent Deng wrote:
>>
>>  > Dear Uwe,
>>  >
>>  > Sorry, I did not describe my question clearly. I created a matrix to
>>  > store color code using rgb function.
>>  >
>>  > abc = rgb(6:36,0,0,maxColorValue = 255)
>>  >
>>  > And after running codes like this
>>  >
>>  > for (i in c(1:20))
>>  > {
>>  >points(...,...,col=abc[i])
>>  >lines(...,col=abc[i])
>>  > }
>>  >
>>  > R still used 8 colors of abc color codes repeatedly to draw the diagram
>>  >
>>  > Any helps?
>>
>>
>>  No, it does not (in fact, all appears to be more or less black on my
>>  screen ;-)). Another example:
>>
>>  plot(1:255, col=rgb(1:255,0,0,maxColorValue = 255))
>>
>>  Uwe Ligges
>>
>>
>>
>>  > Best Regards...
>>  >
>>  > On 12/27/05, Uwe Ligges <[EMAIL PROTECTED]> wrote:
>>  >
>>  >>Vincent Deng wrote:
>>  >>
>>  >>>Hi,
>>  >>>
>>  >>>I'm a new hand in R language. I have about 20 groups of data[x,y] and
>>  >>>want to plot them on a graph. To do this, I write a for-loop as
>>  >>>following: (some codes are omitted for simplicity)
>>  >>>
>>  >>>for (i in c(1:20))
>>  >>>{
>>  >>>  points(...,...,col=i)
>>  >>>  lines(...,col=i)
>>  >>>}
>>  >>>
>>  >>>The problem is "R only plot them with 8 colors repeatly". Could anyone
>>  >>>help me solve this problem? Or is there any package providing plot
>>  >>>function without color limit?
>>  >>
>>  >>
>>  >>After typing
>>  >>
>>  >>  ?colors
>>  >>
>>  >>I get a nice help page that points me to a lot of other functions that
>>  >>generate more than 8 colors. Maybe your installation of R is broken and
>>  >>you cannot see this help page? You certainly tried to get help on colors
>>  >>as well.
>>  >>
>>  >>There is no limit of the color number in the functions above, simply
>>  >>specify the color you want to get. The only color limit applies for the
>>  >>device and for most devices and rgb colors this is 256^3.
>  > >>
>  > >>Uwe Ligges
>  > >>
>  > >>
>  > >>
>  > >>
>  > >>
>  > >>>Best Regards...
>  > >>>
>  > >>>__
>  > >>>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
>>  >>
>>  >>
>>
>>
>
>__
>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


-- 
--
Don MacQueen
Environmental Protection Department
Lawrence Livermore National Laboratory
Livermore, CA, USA

__
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


Re: [R] Help with Kriging

2005-12-29 Thread Paulo Justiniano Ribeiro Jr

Jennifer

The algorithms in geoR are (at least up to now) focused on kriging with 
global neighbourhood.


You you would like to use moving neighbourhood based on distances my nest 
advice is to use a function from another package such as RandomFields, 
gstat, or any other geostats package


Best
P.J.


On Thu, 29 Dec 2005, Jennifer Lenz wrote:


Date: Thu, 29 Dec 2005 11:13:08 -0500
From: Jennifer Lenz <[EMAIL PROTECTED]>
To: r-help@stat.math.ethz.ch
Subject: [R] Help with Kriging

R Experts,

I'm looking for some help with the geoR package.  I'm trying to krig
some data without using a global neighborhood.  I would like to set my
moving neighborhood to a distance, say 100 meters, where I know my data
is spatially correlated.  I have tried the ksline function, but that
only allows my moving neighborhood to be set to a number of data points.
 But, since my data is not equally spaced it makes more sense to use a
distance.  Is there a function out there that would allow me to do this?

Thanks,

Jen Lenz
Tufts University

__
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




Paulo Justiniano Ribeiro Jr
LEG (Laboratório de Estatística e Geoinformação)
Departamento de Estatística
Universidade Federal do Paraná
Caixa Postal 19.081
CEP 81.531-990
Curitiba, PR  -  Brasil
Tel: (+55) 41 3361 3573
Fax: (+55) 41 3361 3141
e-mail: [EMAIL PROTECTED]
http://www.est.ufpr.br/~paulojus__
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

[R] Getting the same y-axis in a multivariate time series plot - plot(ts(...)); ylim does not do the trick

2005-12-29 Thread Søren Højsgaard
I try to obtain the same y-axis for a 2-dim time series with
 
  plot(ts(cbind(rnorm(10), rnorm(10,mean=4))),ylim=c(0,20))

but that does not work. Looking in the code for plot.ts, the ylim-argument 
seems to be taken care of, but not the way I expect. Can anyone help on this?
Thanks
Søren

__
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


Re: [R] search in matrix

2005-12-29 Thread Marc Schwartz (via MN)
On Thu, 2005-12-29 at 17:20 +0100, Florent Bresson wrote:
> I'm dealing with a matrix like :
> 
>  "x"  "y"  "z"
> [1,]  24 1
> [2,]  61 2
> ...
> [n,]  73 1
> 
> For each row I would like to know the header of the
> column which corresponds to the minimum value. In the
> case of my matrix, I would like to obtain the
> following vector :
> 
> z y ... z
> 
> Any idea ?


> mat
 x y z
[1,] 2 4 1
[2,] 6 1 2
[3,] 7 3 1


> colnames(mat)[apply(mat, 1, which.min)]
[1] "z" "y" "z"


apply() (using '1') passes each row to the function which.min() which
returns the index of the minimum value for the row. The indices are then
passed to colnames(mat), returning the column names on a row by row
basis as a vector.

See ?colnames, ?apply and ?which.min for more information.

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


[R] search in matrix

2005-12-29 Thread Florent Bresson
I'm dealing with a matrix like :

 "x"  "y"  "z"
[1,]  24 1
[2,]  61 2
...
[n,]  73 1

For each row I would like to know the header of the
column which corresponds to the minimum value. In the
case of my matrix, I would like to obtain the
following vector :

z y ... z

Any idea ?

__
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


[R] Help with Kriging

2005-12-29 Thread Jennifer Lenz
R Experts,

I'm looking for some help with the geoR package.  I'm trying to krig 
some data without using a global neighborhood.  I would like to set my 
moving neighborhood to a distance, say 100 meters, where I know my data 
is spatially correlated.  I have tried the ksline function, but that 
only allows my moving neighborhood to be set to a number of data points. 
  But, since my data is not equally spaced it makes more sense to use a 
distance.  Is there a function out there that would allow me to do this?

Thanks,

Jen Lenz
Tufts University

__
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


Re: [R] Segmetation Fault in R

2005-12-29 Thread Kort, Eric

Marcelo Damasceno wrote... 

>Hi all,
>
>I has a C code in Linux, it has 7 pointers and compile e run OK, but when I
>run in R happens Segmetation Fault.
>When I use calloc function, it returns NULL.
>What's wrong?
>I would like more information about R-alloc function?
>Thanks!

What is wrong is that there is a bug in your C code related.

A search of the R-help mailing list archives seems to indicate that there have 
been about 10,800 prior posts about segmentation faults, the last of which was 
a few days ago--in which I and then Dr. Ripley described in some detail what 
causes segmentation faults.  I would review these, and if you are still having 
this problem, post a snippet of code that recreates the problem to the 
appropriate C development list serve.

HTH,
Eric

This email message, including any attachments, is for the so...{{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


[R] S4 classes: referencing slots with other slots

2005-12-29 Thread A.J. Rossini
For those who suggest other ways to do this, I ALREADY HAVE ANOTHER
DESIGN SOLUTION, DESCRIBED AT THE END.

That being said, I want to know if it's possible to reference a slot
in an S4 class from another slot, i.e. I'd like to have the "self.*"
semantics of Python so that I can reuse a slot.  That is, for various
reasons it would be nice to be able to do something like:

setClass("fooWfcn",
 representation(dat1="vector",
dat2="vector",
fn1="function",
fn2="function"),
 prototype=list(dat1=0:10,
   dat2=10:20,
   fn1=function(x) { return(x - mean(self.dat1)) },
   fn2=function() { mean(self.dat2) }))

and in the context of

foo <- new("fooWfcn")

have self.dat2 refer to [EMAIL PROTECTED], etc (I.e. in an instantiated object,
have the reference be within the object).

One could easily imagine doing this on the fly, as well, during the "new" call.

 I've looked this up in a number of books, on-line talks/papers, etc,
but havn't managed to find the right page describing it as possible or
impossible.  I think it is the latter (impossible), since one could
easily end up with a nasty self-referencing infinite loop (fn1
refering to fn2 refering to fn1).  If it's the former, I'd be
interested in knowing about it.

Anyway, here's the DESIGN SOLUTION: write methods which do the
appropriate "combination".

Critique:

pro - it works, it's simple, and I've already done it.

con - for the problem I'm looking at, it's not quite so clean, adding
one more layer of indirection that in Python or CLOS I'd not need,
multiplied by a fair number of subclasses.

best,
-tony

[EMAIL PROTECTED]
Muttenz, Switzerland.
"Commit early,commit often, and commit in a repository from which we can easily
roll-back your mistakes" (AJR, 4Jan05).

__
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


Re: [R] Repeating functions

2005-12-29 Thread Sean Davis



On 12/29/05 11:01 AM, "Kort, Eric" <[EMAIL PROTECTED]> wrote:

> Ronnie Babigumira said...
>> 
>> Hi, I have a number of spatial weight files and using Roger Bivand's spdep, I
>> would like to
>> 
>> 1. Convert them into neighbor lists using
>> 2. Convert the neighbor lists into spatial weights
>> 
>> For a given file, the syntax would be
>> 
>> mygal_nb1 <- read.gal("mygalfile1", override.id = TRUE)
>> myweight1 <- nb2listw(mygal_nb1)
>> 
>> I have mygalfile[i] with i from 1 through to 6 and would like to repeat the
>> above two lines through 1 to 6.
>> something like (the syntax below is not correct...just an illustration)
>> 
>> for (i in 1:6) {
>> mygal_nb[i] <- read.gal("mygalfile[i]", override.id = TRUE)
>> myweight[i] <- nb2listw(mygal_nb[i]
>> }
>> 
>> 
>> Kindly point me in the right direction (whilst read through the material I
>> have)
> 
> Looks like you want to create lists of objects.  Likely the best place to
> start is looking at the section on lists in An Introduction to R (which came
> with your R distribution and can be accessed via help.start(), or you can
> download it here: http://cran.r-project.org/doc/manuals/R-intro.pdf).

To be just a bit more explicit, see (untested) code.  Lists are quite
powerful data structures, so learn to use them well.

 mygal_nb <- list()
 myweight <- list()
 for (i in 1:6) {
 mygal_nb[[i]] <- read.gal(mygalfile[i], override.id = TRUE)
 myweight[[i]] <- nb2listw(mygal_nb[[i]])
 }

Sean

__
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


Re: [R] Repeating functions

2005-12-29 Thread Kort, Eric
Ronnie Babigumira said...
>
>Hi, I have a number of spatial weight files and using Roger Bivand's spdep, I 
>would like to
>
>1. Convert them into neighbor lists using
>2. Convert the neighbor lists into spatial weights
>
>For a given file, the syntax would be
>
>mygal_nb1 <- read.gal("mygalfile1", override.id = TRUE)
>myweight1 <- nb2listw(mygal_nb1)
>
>I have mygalfile[i] with i from 1 through to 6 and would like to repeat the 
>above two lines through 1 to 6.
>something like (the syntax below is not correct...just an illustration)
>
>for (i in 1:6) {
>mygal_nb[i] <- read.gal("mygalfile[i]", override.id = TRUE)
>myweight[i] <- nb2listw(mygal_nb[i]
>}
>
>
>Kindly point me in the right direction (whilst read through the material I 
>have)

Looks like you want to create lists of objects.  Likely the best place to start 
is looking at the section on lists in An Introduction to R (which came with 
your R distribution and can be accessed via help.start(), or you can download 
it here: http://cran.r-project.org/doc/manuals/R-intro.pdf).

-Eric

>Ronnie
This email message, including any attachments, is for the so...{{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


Re: [R] Split graph labels in 2 levels

2005-12-29 Thread P Ehlers


Marc Schwartz wrote:

> On Thu, 2005-12-29 at 14:49 +0100, Andrej Kastrin wrote:
> 
>>Dear R users,
>>
>>is there any simple low-level function that split "single-line" graph 
>>labels and produce something like (e.g. for x axis):
>>
>>100300500   700...
>> 200400   600
>>
>>Cheers, Andrej
> 
> 
> You could do something like this:
> 
>  # Draw some points
>  # Do not plot the x axis
>  plot(rnorm(1000), xaxt = "n")
> 
>  # Now create the x axis labels, using "\n" for the odd values
>  # This puts the following even values one line below
>  x.lab <- paste(seq(0, 1000, 100), c("", "\n"), sep = "")
> 
>  # Now do the axis, but tickmarks only
>  axis(1, at = seq(0, 1000, 100), labels = NA)
> 
>  # Now do the labels
>  mtext(1, at = seq(0, 1000, 100), text = x.lab, line = 2)
> 
> 
> See ?axis, ?paste and ?mtext for more information.
> 

Another way is to use the padj argument to axis:

   plot(rnorm(1000), xaxt = "n")
   axis(1, at = seq(0, 1000, 100), padj = c(0, 1))

Adjust padj to suit your needs; e.g. padj = c(1, -0.2)

Peter Ehlers


> You might also want to look at R FAQ 7.27 on rotating axis labels,
> depending upon your requirements.
> 
> 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

__
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


Re: [R] Importing Genstat files into R

2005-12-29 Thread justin bem
Sorry, Graham

I didn't understand first. How cant I get Genstat ?
Best wishes for 2006. 

--- Graham Smith <[EMAIL PROTECTED]> a écrit :

> Justin,
> 
> Thanks for the reply.
> 
> There is no problem getting files out of Genstat and
> into R, indeed Gensat
> exports to R as well as  txt, xls and most other
> statistics program formats.
> It would however it be a lot simpler if R read the
> Genstat files directly.
> 
> I was just rather hoping that someone might have
> written something to do
> this.
> 
> Graham
> 
> On 12/29/05, justin bem <[EMAIL PROTECTED]> wrote:
> >
> > I dont know Genstat, but I if possible to export
> you
> > file to text or other readable format in foriegn
> why
> > not do it ?
> >
> >
> >
> > --- Graham Smith <[EMAIL PROTECTED]> a écrit :
> >
> > > Does anyone know if there is a package or other
> > > method of reading Genstat
> > > files directly into R. Genstat isn't listed in
> the
> > > foreign package.
> > >
> > > Many thanks,
> > >
> > > Graham
> > >
> > >   [[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
> > >
> >
> >
> > Justin BEM
> > Elève Ingénieur Statisticien Economiste
> > BP 294 Yaoundé.
> > Tél (00237)9597295.
> >
> > __
> > 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
> >
>

__
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


Re: [R] Importing Genstat files into R

2005-12-29 Thread Marc Schwartz
Another possibility might be a utility called DataLoad by David Baird.

It is available from:

http://www.american.edu/academic.depts/cas/econ/gaussres/utilitys/dataload.htm

and appears to support the conversion of GenStat files to other formats
(such as CSV), which you could then directly import into R.

HTH,

Marc Schwartz


On Thu, 2005-12-29 at 15:15 +0100, justin bem wrote:
> I dont know Genstat, but I if possible to export you
> file to text or other readable format in foriegn why
> not do it ?
> 
> 
> 
> --- Graham Smith <[EMAIL PROTECTED]> a écrit :
> 
> > Does anyone know if there is a package or other
> > method of reading Genstat
> > files directly into R. Genstat isn't listed in the
> > foreign package.
> > 
> > Many thanks,
> > 
> > Graham

__
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

Re: [R] bVar slot of lmer objects and standard errors

2005-12-29 Thread Doran, Harold
Uli:

The graphic in the paper, sometimes called a catepillar plot, must be
created with some programming as there is (as far as I know) not a
built-in function for such plots. As for the contents of bVar you say
the dimensions are 2,2,28 and there are two random effects and 28
schools. So, from what I know about your model, the third dimension
represents the posterior covariance matrix for each of your 28 schools
as Spencer notes.

For example, consider the following model
> library(Matrix)
> library(mlmRev)
> fm1 <- lmer(math ~ 1 + (year|schoolid), egsingle)

Then, get the posterior means (modes for a GLMM)
> [EMAIL PROTECTED]

These data have 60 schools, so you will see ,,1 through ,,60 and the
elements of each matrix are posterior variances on the diagonals and
covariances in the off-diags (upper triang) corresponding to the
empirical Bayes estimates for each of the 60 schools.

, , 1

   [,1] [,2]
[1,] 0.01007129 -0.001272618
[2,] 0.  0.004588049


Does this help?

Harold


-Original Message-
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] On Behalf Of Spencer Graves
Sent: Thursday, December 29, 2005 6:58 AM
To: Ulrich Keller
Cc: r-help
Subject: Re: [R] bVar slot of lmer objects and standard errors

  Have you received a satisfactory reply to this post?  I
haven't seen one.  Unfortunately, I can't give a definitive answer, but
I can offer an intelligent guess.  With luck, this might encourage
someone who knows more than I do to reply.  If not, I hope these
comments help you clarify the issue further, e.g., by reading the source
or other references.

  I'm not not sure, but I believe that
[EMAIL PROTECTED],,i] is the upper triangular part of the
covariance matrix of the random effects for the i-th level of schoolid. 
  The lower triangle appears as 0, though the code (I believe) iterprets
it as equal to the upper triangle.  More precisely, I suspect it is
created from something that is stored in a more compact form, i.e.,
keeping only a single copy of the off-diagonal elements of symmetric
matrices.  I don't seem to have access to your "nlmframe", so I can't
comment further on those specifics.  You might be able to clarify this
by reading the source code.  I've been sitting on this reply for several
days without finding time to do more with it, so I think I should just
offer what I suspect.

  The specifics of your question suggest to me that you want to
produce something similar to Figure 1.12 in Pinheiro and Bates (2000)
Mixed-Effects Models in S and S-Plus (Springer).  That was produced from
an "lmList" object, not an "lme" object, so we won't expect to get their
exact answers.  Instead, we would hope to get tighter answers available
from pooling information using "lme";  the function "lmList" consideres
each subject separately with no pooling.  With luck, the answers should
be close.

  I started by making a local copy of the data:

library(nlme)
OrthoFem <- Orthodont[Orthodont$Sex=="Female",]

  Next, I believe to switch to "lme4", we need to quit R
completely and restart.  I did that.  Then with the following sequence
of commands I produced something that looked roughly similar to the
confidence intervals produced with Figure 1.12:

library(lme4)
fm1OrthF. <- lmer(distance~age+(age|Subject), data=OrthoFem)

fm1.s <-  coef(fm1OrthF.)$Subject
fm1.s.var <- [EMAIL PROTECTED]
fm1.s0.s <- sqrt(fm1.s.var[1,1,])
fm1.s0.a <- sqrt(fm1.s.var[2,2,])
fm1.s[,1]+outer(fm1.s0.s, c(-2,0,2))
fm1.s[,2]+outer(fm1.s0.a, c(-2,0,2))

  hope this helps.  
  Viel Glueck.
  spencer graves

Ulrich Keller wrote:

> Hello,
> 
> I am looking for a way to obtain standard errors for emprirical Bayes
estimates of a model fitted with lmer (like the ones plotted on page 14
of the document available at
http://www.eric.ed.gov/ERICDocs/data/ericdocs2/content_storage_01/00
0b/80/2b/b3/94.pdf). 


Harold Doran mentioned
(http://tolstoy.newcastle.edu.au/~rking/R/help/05/08/10638.html)
that  the posterior modes' variances can be found in the bVar slot of
lmer objects. However, when I fit e.g. this model:
> 
> lmertest1<-lmer(mathtot~1+(m_escs_c|schoolid),hlmframe)
> 
> then [EMAIL PROTECTED] is a three-dimensional array with
dimensions (2,2,28). 
The factor schoolid has 28 levels, and there are random effects for the
intercept and m_escs_c, but what does the third dimension correspond to?
In other words, what are the contents of bVar, and how can I use them to
get standard errors?
> 
> Thanks in advance for your answers and Merry Christmas,
> 
> Uli Keller
> 
> __
> 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

--
Spencer Graves, PhD
Senior Development Engineer
PDF Solutions, Inc.
333 West San Carlos Street Suite 700
San Jose, CA 95110, USA

[EMAIL PROTECTED]
www.pdf.co

Re: [R] Importing Genstat files into R

2005-12-29 Thread justin bem
I dont know Genstat, but I if possible to export you
file to text or other readable format in foriegn why
not do it ?



--- Graham Smith <[EMAIL PROTECTED]> a écrit :

> Does anyone know if there is a package or other
> method of reading Genstat
> files directly into R. Genstat isn't listed in the
> foreign package.
> 
> Many thanks,
> 
> Graham
> 
>   [[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
> 


Justin BEM
Elève Ingénieur Statisticien Economiste
BP 294 Yaoundé.
Tél (00237)9597295.

__
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


Re: [R] Split graph labels in 2 levels

2005-12-29 Thread Marc Schwartz
On Thu, 2005-12-29 at 14:49 +0100, Andrej Kastrin wrote:
> Dear R users,
> 
> is there any simple low-level function that split "single-line" graph 
> labels and produce something like (e.g. for x axis):
> 
> 100300500   700...
>  200400   600
> 
> Cheers, Andrej

You could do something like this:

 # Draw some points
 # Do not plot the x axis
 plot(rnorm(1000), xaxt = "n")

 # Now create the x axis labels, using "\n" for the odd values
 # This puts the following even values one line below
 x.lab <- paste(seq(0, 1000, 100), c("", "\n"), sep = "")

 # Now do the axis, but tickmarks only
 axis(1, at = seq(0, 1000, 100), labels = NA)

 # Now do the labels
 mtext(1, at = seq(0, 1000, 100), text = x.lab, line = 2)


See ?axis, ?paste and ?mtext for more information.

You might also want to look at R FAQ 7.27 on rotating axis labels,
depending upon your requirements.

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


[R] Repeating functions

2005-12-29 Thread Ronnie Babigumira
Hi, I have a number of spatial weight files and using Roger Bivand's spdep, I 
would like to

1. Convert them into neighbor lists using
2. Convert the neighbor lists into spatial weights

For a given file, the syntax would be

mygal_nb1 <- read.gal("mygalfile1", override.id = TRUE)
myweight1 <- nb2listw(mygal_nb1)

I have mygalfile[i] with i from 1 through to 6 and would like to repeat the 
above two lines through 1 to 6.
something like (the syntax below is not correct...just an illustration)

for (i in 1:6) {
mygal_nb[i] <- read.gal("mygalfile[i]", override.id = TRUE)
myweight[i] <- nb2listw(mygal_nb[i]
}


Kindly point me in the right direction (whilst read through the material I have)

Ronnie

__
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


Re: [R] loop

2005-12-29 Thread Gichangi, Anthony
Your loop will store results of count =6 because
every time the loop executes the results are put in
qw so you replace the previous results.  On the other
hand the results of count=1 is not extracted you basicaly
start at 2.

My advice is you initialize qw before the loop and stack the results
one after the other e.g.

 qw <-NULL
 count<-1
  repeat{
 qw<-c(qw,  c(g[count],1:i[count]) )
  count<-count + 1
  if(count>5) break
 }

Kind regards
Anthony Gichangi

- Original Message - 
From: "gynmeerut" <[EMAIL PROTECTED]>
To: 
Sent: Thursday, December 29, 2005 1:41 PM
Subject: [R] loop


>
> Dear All,
>
> I have to use loop over an array  so I am using following procedure
>
> count<-1
> repeat{
> count<-count + 1
> c(g[count],1:i[count]) ->qw
> if(count>5)break
> }
>
> as  a result qw is
> [1]  0.9643836  1.000  2.000  3.000  4.000  5.000
> [7]  6.000  7.000  8.000  9.000 10.000 11.000
> [13] 12.000 13.000 14.000 15.000 16.000 17.000
> [19] 18.000 19.000 20.000 21.000 22.000 23.000
> [25] 24.000 25.000 26.000 27.000 28.000 29.000
> [31] 30.000 31.000 32.000 33.000 34.000 35.000
> [37] 36.000 37.000 38.000 39.000
>
> which is according to the last value of  count i.e 6. What is the problem 
> in error
>
>
> but I expect the result as
>
> [1]  0.8328767  0.2410959  0.5315068  0.1424658  0.2520548  0.9643836
> [7]  6.000  7.000  8.000  0.2410959  1.000  2.000
> [13]  3.000  4.000  5.000  0.5315068  1.000  2.000
> [19]  3.000  4.000  5.000  6.000  7.000  8.000
> [25]  9.000 10.000 11.000 12.000 13.000  0.1424658
> [31]  1.000  0.000  0.2520548  1.000  2.000  3.000
> [37]  4.000  5.000  6.000  7.000  8.000  9.000
> [43] 10.000 11.000 12.000 13.000 14.000 15.000
> [49] 16.000 17.000 18.000 19.000 20.000 21.000
> [55] 22.000 23.000  0.9643836  1.000  2.000  3.000
> [61]  4.000  5.000  6.000  7.000  8.000  9.000
> [67] 10.000 11.000 12.000 13.000 14.000 15.000
> [73] 16.000 17.000 18.000 19.000 20.000 21.000
> [79] 22.000 23.000 24.000 25.000 26.000 27.000
> [85] 28.000 29.000 30.000 31.000 32.000 33.000
> [91] 34.000 35.000 36.000 37.000 38.000 39.000
>
> which is combination of 6 vectors of lengths {9,6,14,3,24,40).
>
>
> Please tell me the error or alternate procedure for looping.
>
> regards,
>
> GS
>
> __
> 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

__
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


[R] Split graph labels in 2 levels

2005-12-29 Thread Andrej Kastrin
Dear R users,

is there any simple low-level function that split "single-line" graph 
labels and produce something like (e.g. for x axis):

100300500   700...
 200400   600

Cheers, Andrej

__
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


[R] Segmetation Fault in R

2005-12-29 Thread Marcelo Damasceno
Hi all,

I has a C code in Linux, it has 7 pointers and compile e run OK, but when I
run in R happens Segmetation Fault.
When I use calloc function, it returns NULL.
What's wrong?
I would like more information about R-alloc function?
Thanks!

--
Marcelo Damasceno de Melo
Graduando em Ciência da Computação
Departamento de Tecnologia da Informação - TCI
Universidade Federal de Alagoas - UFAL
Maceió - Alagoas - Brasil
Projeto CoCADa - Construção do Conhecimento por Agrupamento de dados
+55 82 8801-2119

[[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

Re: [R] loop

2005-12-29 Thread Duncan Murdoch
On 12/29/2005 7:41 AM, gynmeerut wrote:
> Dear All,
> 
> I have to use loop over an array  so I am using following procedure
>  
> count<-1
>  repeat{
>  count<-count + 1
>  c(g[count],1:i[count]) ->qw
>  if(count>5)break
>  }

We can't reproduce this, as we don't have g or i.  But the general 
advice in a case like this is to simulate the loop by hand:  write down 
what is in each of the variables, and walk through the loop.

I suspect your problem is in initializing count improperly, or in 
putting the test in the wrong place.

Duncan Murdoch
> 
> as  a result qw is
> [1]  0.9643836  1.000  2.000  3.000  4.000  5.000
>  [7]  6.000  7.000  8.000  9.000 10.000 11.000
> [13] 12.000 13.000 14.000 15.000 16.000 17.000
> [19] 18.000 19.000 20.000 21.000 22.000 23.000
> [25] 24.000 25.000 26.000 27.000 28.000 29.000
> [31] 30.000 31.000 32.000 33.000 34.000 35.000
> [37] 36.000 37.000 38.000 39.000
> 
> which is according to the last value of  count i.e 6. What is the problem in 
> error
>  
> 
> but I expect the result as
> 
>  [1]  0.8328767  0.2410959  0.5315068  0.1424658  0.2520548  0.9643836
>  [7]  6.000  7.000  8.000  0.2410959  1.000  2.000
> [13]  3.000  4.000  5.000  0.5315068  1.000  2.000
> [19]  3.000  4.000  5.000  6.000  7.000  8.000
> [25]  9.000 10.000 11.000 12.000 13.000  0.1424658
> [31]  1.000  0.000  0.2520548  1.000  2.000  3.000
> [37]  4.000  5.000  6.000  7.000  8.000  9.000
> [43] 10.000 11.000 12.000 13.000 14.000 15.000
> [49] 16.000 17.000 18.000 19.000 20.000 21.000
> [55] 22.000 23.000  0.9643836  1.000  2.000  3.000
> [61]  4.000  5.000  6.000  7.000  8.000  9.000
> [67] 10.000 11.000 12.000 13.000 14.000 15.000
> [73] 16.000 17.000 18.000 19.000 20.000 21.000
> [79] 22.000 23.000 24.000 25.000 26.000 27.000
> [85] 28.000 29.000 30.000 31.000 32.000 33.000
> [91] 34.000 35.000 36.000 37.000 38.000 39.000
> 
> which is combination of 6 vectors of lengths {9,6,14,3,24,40).
> 
> 
> Please tell me the error or alternate procedure for looping.
> 
> regards,
> 
> GS
> 
> __
> 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

__
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


Re: [R] Open a new script from R command prompt

2005-12-29 Thread Ronnie Babigumira
Hi Justin
Thanks. Kort had pointed me in the same direction however, using this means I 
cannot see my output until I close the 
editor window. Bodgan proposed Alt-F-N and while this is not a command issued 
at the command prompt, it opens the editor 
without having to use the mouse so I will live with it.

Ronnie



justin bem wrote:
>   Try edit() or ?edit
> 
> Ronnie Babigumira <[EMAIL PROTECTED]> a écrit :   Hi, (this is a minor 
> irritation), is it possible for me to call R's editor from the R command 
> prompt (I searched for 
> script but that didn't yield anything). (The mouse-file-new-script route is a 
> minor irritation )
> 
> __
> 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
> 
> 
> 
> 
>   
> -
> 
>   [[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

__
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


[R] loop

2005-12-29 Thread gynmeerut

Dear All,

I have to use loop over an array  so I am using following procedure
 
count<-1
 repeat{
 count<-count + 1
 c(g[count],1:i[count]) ->qw
 if(count>5)break
 }

as  a result qw is
[1]  0.9643836  1.000  2.000  3.000  4.000  5.000
 [7]  6.000  7.000  8.000  9.000 10.000 11.000
[13] 12.000 13.000 14.000 15.000 16.000 17.000
[19] 18.000 19.000 20.000 21.000 22.000 23.000
[25] 24.000 25.000 26.000 27.000 28.000 29.000
[31] 30.000 31.000 32.000 33.000 34.000 35.000
[37] 36.000 37.000 38.000 39.000

which is according to the last value of  count i.e 6. What is the problem in 
error
 

but I expect the result as

 [1]  0.8328767  0.2410959  0.5315068  0.1424658  0.2520548  0.9643836
 [7]  6.000  7.000  8.000  0.2410959  1.000  2.000
[13]  3.000  4.000  5.000  0.5315068  1.000  2.000
[19]  3.000  4.000  5.000  6.000  7.000  8.000
[25]  9.000 10.000 11.000 12.000 13.000  0.1424658
[31]  1.000  0.000  0.2520548  1.000  2.000  3.000
[37]  4.000  5.000  6.000  7.000  8.000  9.000
[43] 10.000 11.000 12.000 13.000 14.000 15.000
[49] 16.000 17.000 18.000 19.000 20.000 21.000
[55] 22.000 23.000  0.9643836  1.000  2.000  3.000
[61]  4.000  5.000  6.000  7.000  8.000  9.000
[67] 10.000 11.000 12.000 13.000 14.000 15.000
[73] 16.000 17.000 18.000 19.000 20.000 21.000
[79] 22.000 23.000 24.000 25.000 26.000 27.000
[85] 28.000 29.000 30.000 31.000 32.000 33.000
[91] 34.000 35.000 36.000 37.000 38.000 39.000

which is combination of 6 vectors of lengths {9,6,14,3,24,40).


Please tell me the error or alternate procedure for looping.

regards,

GS

__
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


[R] calculating recursive sequences

2005-12-29 Thread Vaidotas Zemlys
Hi,


I was trying to repeat the estimation of threshold GARCH models from
the book "Analysis of Financial Time Series" by Ruey S. Tsay, and I
was succesfull, but I had to use "for" loop, which is quite slow. The
loop is necessary, since you need to calculate recursive sequence. Is
there a faster way to do this in R, without using loops?

The model is such:
r_t = \mu + \alpha_2 r_{t-2} + a_t
a_t = \sigma_t\varepsilon_t

\sigma_t^2 =
\beta_1a_{t-1}^2+\beta_2\sigma_{t-1}^2+
1_{\{a_{t-1}>0\}}(\gamma_0+
\gamma_1a_{t-1}^2+\gamma_2\sigma^2_{t-1})

It is asummed that \varepsilon_t are iid and normal with zero mean and
variance one. The data given is r_t, and you have to estimate
variables, \mu, \alpha, \beta and \gamma. Since

\varepsilon_t=\frac{a_t}{\sqrt{sigma_t}}

using the equations we calculate a_t and \sigma_t and estimate the
variables using maximum likelihood method. a_t can be estimated
directly using first equation and rt. \sigma_t^2 depends on
sigma_{t-1}^2, so it should be calculated recursively.

The function calculating negative log-likelihood of this problem I wrote:

garchln <- function(p,rt) {
n <- length(rt)

at <- rt[4:n]-p[1]-p[2]*rt[4:n-2]
u <- as.numeric(at>0)
h <- rep(0,length(at))
  # h is \sigma_t^2
   for(i in 1:(length(h)-1)) {
h[i+1] <- p[3]*at[i]^2+p[4]*h[i]+u[i]*(p[5]+p[6]*at[i]^2+p[7]*h[i])
}

   #Maximum likelihood function
sum(log(h[-1])+(at[-1]^2)/h[-1])/2
#list(h=h[-1],at=at[-1])
}

For fitting I used optim, with methods "Nelder-Mead" and "BFGS",

Initial parameter values from the book are
 0.03 -0.03  0.10  0.60  0.10  0.05  0.10
The fitted values from the book are
 0.043 -0.022  0.098  0.954  0.060 -0.052 -0.069.

The link to the data used:
http://www.gsb.uchicago.edu/fac/ruey.tsay/teaching/fts/d-ibmln99.dat

For this problem recursive sequence is linear, so it is possible to
calculate it as a linear equations solution, but it is easy to think
of the case where the recursion is non-linear. Is the  speed-up
possible only by writing C or Fortran code with loops?

Vaidotas Zemlys
--
Doctorate student, http://www.mif.vu.lt/katedros/eka/katedra/zemlys.php
Vilnius University

__
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


Re: [R] Testing a linear hypothesis after maximum likelihood

2005-12-29 Thread Spencer Graves
  Why can't you use a likelihood ratio?  I would write two slightly 
different functions, the second of which would use the linear constraint 
to eliminate one of the coefficients.  Then I'd refer 2*log(likelihood 
ratio) to chi-square(1).  If I had some question about the chi-square 
approximation to the distribution of that 2*log(likelihood ratio) 
statistic, I'm use some kind of Monte Carlo, e.g., MCMC.

  If you'd like more help from this listserve, PLEASE do read the 
posting guide! "www.R-project.org/posting-guide.html".  Anecdotal 
evidence suggests that posts that follow more closely the suggestions in 
that guide tend to get more useful replies quicker.

  hope this helps.
  spencer graves


Peter Muhlberger wrote:

> I'd like to be able to test linear hypotheses after setting up and running a
> model using optim or perhaps nlm.  One hypothesis I need to test are that
> the average of several coefficients is less than zero, so I don't believe I
> can use the likelihood ratio test.
> 
> I can't seem to find a provision anywhere for testing linear combinations of
> coefficients after max. likelihood.
> 
> Cheers & happy holidays,
> 
> Peter
> 
> __
> 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

-- 
Spencer Graves, PhD
Senior Development Engineer
PDF Solutions, Inc.
333 West San Carlos Street Suite 700
San Jose, CA 95110, USA

[EMAIL PROTECTED]
www.pdf.com 
Tel:  408-938-4420
Fax: 408-280-7915

__
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


Re: [R] bVar slot of lmer objects and standard errors

2005-12-29 Thread Spencer Graves
  Have you received a satisfactory reply to this post?  I haven't seen 
one.  Unfortunately, I can't give a definitive answer, but I can offer 
an intelligent guess.  With luck, this might encourage someone who knows 
more than I do to reply.  If not, I hope these comments help you clarify 
the issue further, e.g., by reading the source or other references.

  I'm not not sure, but I believe that
[EMAIL PROTECTED],,i] is the upper triangular part of the 
covariance matrix of the random effects for the i-th level of schoolid. 
  The lower triangle appears as 0, though the code (I believe) iterprets 
it as equal to the upper triangle.  More precisely, I suspect it is 
created from something that is stored in a more compact form, i.e., 
keeping only a single copy of the off-diagonal elements of symmetric 
matrices.  I don't seem to have access to your "nlmframe", so I can't 
comment further on those specifics.  You might be able to clarify this 
by reading the source code.  I've been sitting on this reply for several 
days without finding time to do more with it, so I think I should just 
offer what I suspect.

  The specifics of your question suggest to me that you want to produce 
something similar to Figure 1.12 in Pinheiro and Bates (2000) 
Mixed-Effects Models in S and S-Plus (Springer).  That was produced from 
an "lmList" object, not an "lme" object, so we won't expect to get their 
exact answers.  Instead, we would hope to get tighter answers available 
from pooling information using "lme";  the function "lmList" consideres 
each subject separately with no pooling.  With luck, the answers should 
be close.

  I started by making a local copy of the data:

library(nlme)
OrthoFem <- Orthodont[Orthodont$Sex=="Female",]

  Next, I believe to switch to "lme4", we need to quit R
completely and restart.  I did that.  Then with the following sequence 
of commands I produced something that looked roughly similar to the 
confidence intervals produced with Figure 1.12:

library(lme4)
fm1OrthF. <- lmer(distance~age+(age|Subject), data=OrthoFem)

fm1.s <-  coef(fm1OrthF.)$Subject
fm1.s.var <- [EMAIL PROTECTED]
fm1.s0.s <- sqrt(fm1.s.var[1,1,])
fm1.s0.a <- sqrt(fm1.s.var[2,2,])
fm1.s[,1]+outer(fm1.s0.s, c(-2,0,2))
fm1.s[,2]+outer(fm1.s0.a, c(-2,0,2))

  hope this helps.  
  Viel Glueck.
  spencer graves

Ulrich Keller wrote:

> Hello,
> 
> I am looking for a way to obtain standard errors for emprirical Bayes 
estimates of a model fitted with lmer (like the ones plotted on page 14
of the document available at
http://www.eric.ed.gov/ERICDocs/data/ericdocs2/content_storage_01/000b/80/2b/b3/94.pdf).
 


Harold Doran mentioned
(http://tolstoy.newcastle.edu.au/~rking/R/help/05/08/10638.html)
that  the posterior modes' variances can be found in the bVar slot of lmer
objects. However, when I fit e.g. this model:
> 
> lmertest1<-lmer(mathtot~1+(m_escs_c|schoolid),hlmframe)
> 
> then [EMAIL PROTECTED] is a three-dimensional array with dimensions (2,2,28). 
The factor schoolid has 28 levels, and there are random effects for the
intercept and m_escs_c, but what does the third dimension correspond to?
In other words, what are the contents of bVar, and how can I use them to
get standard errors?
> 
> Thanks in advance for your answers and Merry Christmas,
> 
> Uli Keller
> 
> __
> 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

-- 
Spencer Graves, PhD
Senior Development Engineer
PDF Solutions, Inc.
333 West San Carlos Street Suite 700
San Jose, CA 95110, USA

[EMAIL PROTECTED]
www.pdf.com 
Tel:  408-938-4420
Fax: 408-280-7915

__
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


Re: [R] Open a new script from R command prompt

2005-12-29 Thread justin bem
  Try edit() or ?edit

Ronnie Babigumira <[EMAIL PROTECTED]> a écrit :   Hi, (this is a minor 
irritation), is it possible for me to call R's editor from the R command prompt 
(I searched for 
script but that didn't yield anything). (The mouse-file-new-script route is a 
minor irritation )

__
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





-

[[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

Re: [R] Open a new script from R command prompt

2005-12-29 Thread justin bem


Ronnie Babigumira <[EMAIL PROTECTED]> a écrit :  Hi, (this is a minor 
irritation), is it possible for me to call R's editor from the R command prompt 
(I searched for 
script but that didn't yield anything). (The mouse-file-new-script route is a 
minor irritation )

__
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




-

[[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

Re: [R] segmetation fault and abort trap

2005-12-29 Thread Prof Brian Ripley
These are C programming questions, and the messages are specific to your 
OS.  Please use a more appropriate list: the posting guide does say
`questions involving C' should go elsewhere.

For the specific question here, I think you need to ask on a MacOS list.

BTW, Eric's definition is a good part of the story but not always all the 
story.  Some systems differentiate bus errors from segfaults, and some do 
not.  Generally segfaults come from SIGSEGV signals, bus errors from 
SIGBUS signals and abort traps from SIGABRT signals (and usually from 
the run-time support code calling the C function abort()).  But I get the 
feeling that is not the level of explanation you are seeking.

On Wed, 28 Dec 2005, Elizabeth Lawson wrote:

> Sorry, one more thing.  Sometimes if crashes with Abort trap and sometimes it 
> crashes with Segmentation fault.  Help!
>
>  Thanks!
>
> Elizabeth Lawson <[EMAIL PROTECTED]> wrote:
>  I fixed the problem that I was having with the segmentation error and can 
> load it with dyn.load.
>
> I have one more problem now that I don't understand. I am using the .C in a 
> for loop
> for( i in 1:ncol(L.D)){
> new[,i]<-.C("mycode",as.double(L.D[,i]))
> }
>
> It crashes at different places each time I try to run it. So I tried
> i<-1
> then running
> new[,i]<-.C("mycode",as.double(L.D[,i]))
> repeatedly.
>
> It works for 3-6 times then will crash with the error
> abort trap
>
> What does this mean?
> Why does it work sometimes and not others?
>
> Any suggestions?
>
> Thanks,
>
> Elizabeth Lawson
> "Kort, Eric" wrote:
>
> Elizabeth Lawson Wrote:
>>
>> Hey,
>>
>> I don;t know if anyone has come across this error before...
>
> More times than I care to remember.
>
>>
>> I am running R on the terminal of my MAC OS X 10.3.4 and I have
> written
>> C code and compiled it using
>> R CMD SHLIB mycode.c
>> There were no problems in compiling so I now have mycode.o and
> mycode.so.
>
> A segmentation fault occurs when you try to access memory you didn't
> allocate for your use (see below). In other words, you are trying to
> use memory outside the segment allocated to your program (which, if
> allowed, could result in corrupting memory being used by other
> programs--which brings back unhappy memories of days gone by when an
> error in one program could crash the whole system). So these kinds of
> problems will not show up at compile time...only when you actually run
> the program.
>
>> I used dyn.load("mycode.so") and again, no problems. But when I try
> to
>> use the code .C("mycode",x) I get the error Segmentation fault and R
> shuts
>> down. Any ideas as to where my problem could be?
>
> When I run into a segmentation fault, it is usually because I am trying
> to access an element of an array that is beyond what I have allocated,
> as in...
>
> int main()
> {
> int *a;
> a = (int*) malloc(3*sizeof(int));
> printf("Fasten your seatbelts...\n");
> a[4000] = 12;
> return(0);
> }
>
> One less obvious way this can happen is forgetting to initialize your
> arrays to the proper length before passing a reference to them to your C
> function.
>
> If you still are having trouble, you could post a small snippet of code
> that recreates the error for us to examine.
>
> HTH,
> Eric
>
>
> This email message, including any attachments, is for the so...{{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
>
>
>
>
> -
>
>
>   [[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
>

-- 
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


Re: [R] Show graph integrated to GUI

2005-12-29 Thread Prof Brian Ripley
This can be done with tkrplot package, whose sole example is of a slider 
and a graph in the same window.

(AFAICS the slider function in the TeachingDemos package shows no more 
concepts than the tkdensity demo which ships with R.)

On Thu, 29 Dec 2005, Chihiro Kuraya wrote:

> Hi,
>
> Thank you for information..
>
> I tried the TechingDemos package.
> But it seems that the window which contains slider control
> and graph window is separated.
>
> What I want to do is show slider control and graph
> in one window simultaneiously.
> It is possible ?
>
> Chihiro Kuraya
>
>
> "Gregory Snow" <[EMAIL PROTECTED]> wrote:
>
>> The slider function in the TeachingDemos and relax packages (same
>> function is in both packages you can use either) provides a way to do
>> this using a Tk window.  There are also several functions in the
>> TechingDemos package that use a lower level interface to a Tk window
>> that you can look at their source code as an example to build your own
>> (examples include: vis.gamma, rotate.persp, and run.power.examp).
>>
>> Hope this helps,

>>> -Original Message-
>>> From: [EMAIL PROTECTED]
>>> [mailto:[EMAIL PROTECTED] On Behalf Of Chihiro Kuraya
>>> Sent: Saturday, December 24, 2005 9:53 PM
>>> To: r-help@stat.math.ethz.ch
>>> Subject: [R] Show graph integrated to GUI
>>>
>>> Hi all,
>>>
>>> It it posssible to show graph which is integrated to some GUI
>>> (e.g. TclTk or R-wxPython).
>>>
>>> I want to make an application by R,
>>> for example, like the following picture:
>>> http://www.natch.co.uk/downloads/SigJenny/SJnScreenShot.gif

-- 
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