[R] apply(ing) to sum subset of a vector

2006-03-27 Thread Fred J.
Dear R users
 
 I am trying to sum selective elements of a vector but my solution
 is not cutting it.
 
 Example:
  g - 1:5;
 
  from - 1:3;
  to - 3:5;
 from to
 1   3
 2   4
 3   5
 
 so I expect 3 sums from g
 1+2+3  that is 1 to 3 of g
 2+3+4  that is 2 to 4 of g
 3+4+5  that is 3 to 5 of g
 
 my solution will not work.
 sum.em - function(g, c1, c2) sum(g[c1:c2])
 apply(g, 1, sum.em, ...) I don't think so because apply is not
 aware of the from and to. and if I f - list(g, from, to) that will not fit 
with the second arg of apply.
 
 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


[R] seq(2,5,-2) not an error but NULL

2006-03-27 Thread Christian Hoffmann
Hi,

This may belong more to r-develop, but general discussion may be useful 
(for the how many-th time ?)

seq(2,5,-2)
seq(5,2,2)

both result in

Error in seq.default(2, 5, -2) : wrong sign in 'by' argument

But often, if not always, mathematicians and programmers want a 
behaviour e.g. in for loops, where this statement results in an empty 
statement, that is

for (ii in seq(2,5,-2)) print(ii)

were equivalent to

for (ii in NULL) print(ii).

The relevant part in seq.default is now

 if (n  0)
 stop(wrong sign in 'by' argument)

but could be changed by option to

   return(NULL)

I think there should be an option to seq requiring this behaviour, or a 
specific function, may be even a special operator, e.g. %;%:

3;5 resulting in NULL.

What do you think?

Christian
-- 
Dr. Christian W. Hoffmann,
Swiss Federal Research Institute WSL
Mathematics + Statistical Computing
Zuercherstrasse 111
CH-8903 Birmensdorf, Switzerland

Tel +41-44-7392-277  (office)   -111(exchange)
Fax +41-44-7392-215  (fax)
[EMAIL PROTECTED]
http://www.wsl.ch/staff/christian.hoffmann

International Conference 5.-7.6.2006 Ekaterinburg Russia
Climate changes and their impact on boreal and temperate forests
http://ecoinf.uran.ru/conference/

__
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] apply(ing) to sum subset of a vector

2006-03-27 Thread jim holtman
create a matrix and then use apply:

 g - 1:5;

  from - 1:3;
  to - 3:5;

 index - cbind(from,to)
 apply(index, 1, function(x) sum(g[x[1]:x[2]]))
[1]  6  9 12




On 3/27/06, Fred J. [EMAIL PROTECTED] wrote:

 Dear R users

 I am trying to sum selective elements of a vector but my solution
 is not cutting it.

 Example:
  g - 1:5;

  from - 1:3;
  to - 3:5;
 from to
 1   3
 2   4
 3   5

 so I expect 3 sums from g
 1+2+3  that is 1 to 3 of g
 2+3+4  that is 2 to 4 of g
 3+4+5  that is 3 to 5 of g

 my solution will not work.
 sum.em - function(g, c1, c2) sum(g[c1:c2])
 apply(g, 1, sum.em, ...) I don't think so because apply is not
 aware of the from and to. and if I f - list(g, from, to) that will not
 fit with the second arg of apply.

 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




--
Jim Holtman
Cincinnati, OH
+1 513 646 9390 (Cell)
+1 513 247 0281 (Home)

What the problem you are trying to solve?

[[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] load huge image

2006-03-27 Thread Henrik Bengtsson
On 3/27/06, Martin Maechler [EMAIL PROTECTED] wrote:
  Gottfried == Gottfried Gruber [EMAIL PROTECTED]
  on Sun, 26 Mar 2006 10:27:35 +0200 writes:

 Gottfried hello, i have run around 65000 regressions and
 Gottfried stored them in a list. then i stored the session
 Gottfried with save.image on my hard disk. the file is
 Gottfried almost 1GB. when i now want to load the image it
 Gottfried took tons of time. even after 12h of loading it
 Gottfried was not done, although the saving was done fairly
 Gottfried fast.

 I'm sure it takes so lang because you (i.e. R) run out of RAM
 and the machine starts to swap.

 Try to get access to a Linux (or other Unix-alike) machine with
 a 64-bit version of R and about 8 GB of RAM (maybe 4 GB is
 already sufficient). I guess then you should be able to read it
 much more quickly.

 For 65000 regressions, do you need more than the estimated
 coefficients or -- a bit more informatively -- the

   coef(summary( lm-fit ))  result ?

 If you had only saved these coefficient matrices, I'm sure you'd
 have need **much** less memory.


 Gottfried i fear i have to run the regressions again and
 Gottfried store them in a database ...

 or really store what you need instead of everything ...

 Gottfried can i load this file? any suggestions?

Do you need to have all of them in memory at once?  Instead of using
save.image() can't you use save()/load() on each of the regression
fits?  You can name the files using sprintf(regression%05d.Rdata,
idx) or similar.

Also, as Martin says, fit objects contains a lot of information that
you might not need; remove these before saving by setting the elements
you don't want to NULL.

/Henrik

 Gottfried thanks  bets regards, gg --
 Gottfried ---
 Gottfried Gottfried Gruber
 Gottfried mailto:[EMAIL PROTECTED] www:
 Gottfried http://gogo.sehrsupa.net

 Gottfried __
 Gottfried R-help@stat.math.ethz.ch mailing list
 Gottfried https://stat.ethz.ch/mailman/listinfo/r-help
 Gottfried PLEASE do read the posting guide!
 Gottfried 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




--
Henrik Bengtsson
Mobile: +46 708 909208 (+2h UTC)

__
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] Missing Argument in optim()

2006-03-27 Thread voodooochild
Hello everybody,

i already searched the archieves, but i still don't know what is wrong 
in my implementation, mybe anybody coud give me some advice

ll1-function(rho,theta,beta1,beta2,beta3,beta4,t,Szenariosw5,Testfaellew5,X1,X2)
 
{
n-length(t) 
t-cumsum(t)
  tn-t[length(t)]
  Szenn-Szenariosw5[length(Szenariosw5)]
  Testn-Testfaellew5[length(Testfaellew5)]
  X1n-X1[length(X1)]
X2n-X2[length(X2)]
n/rho-(tn^theta)*(beta1*Szenn+beta2*Testn+beta3*X1n+beta4*X2n)
}

p0-c(rho=1,theta=1,beta1=1,beta2=1,beta3=1,beta4=1)

optim(p0,fn=ll1,t=t,Szenariosw5=Szenariosw5,Testfaellew5=Testfaellew5,X1=X1,X2=X2)

i always got an error message: Argument theta is missing, with no default
but i don't know what is wrong?

thanks!

__
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 create a directoy with R

2006-03-27 Thread pau carre
Hello, I am trying to create directories with R. I would like R to
create directories because it is platform independent. I tried using
file() and searching in R Data Import/Export but I did not succeed.
I think it must be some function since exists the unlink to remove
directories (and files).

Pau

__
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] Missing Argument in optim()

2006-03-27 Thread Prof Brian Ripley
On Mon, 27 Mar 2006, [EMAIL PROTECTED] wrote:

 Hello everybody,

 i already searched the archieves, but i still don't know what is wrong
 in my implementation, mybe anybody coud give me some advice

 ll1-function(rho,theta,beta1,beta2,beta3,beta4,t,Szenariosw5,Testfaellew5,X1,X2)
 {
n-length(t)
t-cumsum(t)
  tn-t[length(t)]
  Szenn-Szenariosw5[length(Szenariosw5)]
  Testn-Testfaellew5[length(Testfaellew5)]
  X1n-X1[length(X1)]
X2n-X2[length(X2)]
n/rho-(tn^theta)*(beta1*Szenn+beta2*Testn+beta3*X1n+beta4*X2n)
 }

 p0-c(rho=1,theta=1,beta1=1,beta2=1,beta3=1,beta4=1)

 optim(p0,fn=ll1,t=t,Szenariosw5=Szenariosw5,Testfaellew5=Testfaellew5,X1=X1,X2=X2)

 i always got an error message: Argument theta is missing, with no default
 but i don't know what is wrong?

From the help page:

   fn: A function to be minimized (or maximized), with first
   argument the vector of parameters over which minimization is
   to take place.  It should return a scalar result.

ll1 is not such a function, and I think you need to replace the first 6 
arguments by a single 6-element vector.

-- 
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] Clustering question \ dist(datmat)

2006-03-27 Thread Sean Davis



On 3/27/06 12:19 AM, kumar zaman [EMAIL PROTECTED] wrote:

 Dear Gabor and all ;

   I know this will work; but i already have a distance matrix calculated using
 my distance measure Dij = 0.5 * ( 1 - cos(theta_i - theta_j)), if i do
 hclust(as.dist(df)) then i am taking distance another time for a matrix  df 
 which is supposed to be a distance matrix, i hope i am clear ;

   ps: I just found out i can use  kmeans(df, 3, iter.max=100) it will take
 df as calculated by Dij. I still need to use methods in hclust like  single,
 average, ward, median, mcquitty, ...etc

   Thank u anyway.

Kumar,

If I understand Your point, you are misunderstanding what as.dist() does.
It does not compute a distance matrix.  Instead, it simply makes a matrix
into a dist object, which is NOT just a matrix.  However, the distances in
a matrix converted to a dist object are not altered.  Therefore, you are
not taking distance another time; instead, you are simply converting the
distance matrix into a form that hclust can understand.

Hope that helps clarify.

Sean


 
 Gabor Grothendieck [EMAIL PROTECTED] wrote:
   A distance matrix must be of class dist. Try
 
 hclust(as.dist(df))
 
 
 On 3/26/06, kumar zaman wrote:
 Hello everybody. I am trying to cluster circular data (data points which are
 angles), thus i can not use the dist function in mclust to generate my
 distance matrix, I am using the function  Dij = 0.5*( 1 - cos(theta_i -
 theta_j)). The thing is hclust will not accept this distance matrix, i
 tried to put it in a data frame, but again i get an error message saying 
 Error in if (n  2) stop(must have n = 2 objects to cluster) : argument is
 of length zero. The distance matrix dist producing is a lower triangular
 one, mine is a square matrix, which i think does not matter. My question how
 to make hclust process my distance matrix, what i am doing wrong. I am sure
 the problem is with the distance matrix format, Any suggestions are highly
 apprciated, the code below shows what i have done.
 
 clust1- as.vector(rvm(5,5,15))
 clust2- as.vector(rvm(5,10,15))
 clust3- as.vector(rvm(5,15,15))
 clust4- as.vector(rvm(5,20,15))
 clust5- as.vector(rvm(5,25,15))
 data1- rbind(clust1,clust2,clust3,clust4,clust5)
 datmat- matrix(data1,nrow=25,ncol=1,byrow=TRUE)
 circ.plot(datmat)
 df- array(dim=c(25,25))
 for (i in 1:25){
 for (j in 1:25){
 df[i,j]- 0.5*(1 - cos(datmat[i] - datmat[j]))
 }
 }
 hcA-hclust(df,method=average)
 
 Ahmed
 Florida
 
 
 -
 
 [[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
 
 
 
 
 Ahmed Albatineh,PhD
 Assistant Professor of Statistics
 Nova Southeastern University
 Fort Lauderdale, FL 33314
 U.S.A
 
 -
 
 [[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] How to create a directoy with R

2006-03-27 Thread David Whiting
?dir.create

On Mon, 2006-03-27 at 13:07 +0200, pau carre wrote:
 Hello, I am trying to create directories with R. I would like R to
 create directories because it is platform independent. I tried using
 file() and searching in R Data Import/Export but I did not succeed.
 I think it must be some function since exists the unlink to remove
 directories (and files).
 
 Pau
 
 __
 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] A plotting question - how to get error bars?

2006-03-27 Thread toby
Dear R list,

Can anyone help with a plotting question? I'm trying to display some data
on a plot and I've almost got the format I need (see code below), but 2
things I can't get:

1. How to get Jan,Feb,Mar on the x=axis instead of 1:3?
2. How to get Ts on the end of my error bars like you have in standard
scientific plots?

Any comments gratefully received!

Thanks,
Toby

xvals=1:3   #couldn't get it to be Jan, Feb, Mar on the x-axis
rgrT1=c(10,20,30)
errbarsT1lo=c(0.5830952,0.3741657,0.8944272)
errbarsT1up=errbarsT1lo
rgrT2=c(25,30,35)
errbarsT2lo=c(1.356466,3.535534,1.140175)
errbarsT2up=errbarsT2lo
minx=min(xvals);maxx=max(xvals)
miny=min(rgrT1-errbarsT1lo,rgrT2-errbarsT2lo);maxy=max(rgrT1+errbarsT1up,rgrT2+errbarsT2up)
plot(x=0,y=0,type=n,xlim=c(minx,maxx),ylim=c(miny,maxy),lab=c(2,20,0),bty=l,xlab=month,ylab=Relative
Growth Rate)
points(x=xvals,y=rgrT1,pch=21)
symbols(x=xvals,y=rgrT1,boxplots=cbind(0,0,errbarsT1lo,errbarsT1up,0.5),inches=FALSE,add=TRUE)
#symbols does the error bars, but without the Ts at the end. 
The
boxplot command does the Ts, but you can't have them without the box in
the middle (and you can't have different symbols for points either)
lines(x=xvals,y=rgrT1,lty=21)
points(x=xvals,y=rgrT2,pch=24)
symbols(x=xvals,y=rgrT2,boxplots=cbind(0,0,errbarsT2lo,errbarsT2up,0.5),inches=FALSE,add=TRUE)
lines(x=xvals,y=rgrT2,lty=24)
legend(x=right,c(Treatment 1,Treatment 2),pch=c(21,24))

__
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 create a directoy with R

2006-03-27 Thread Sarah Goslee
I think you need to use system(mkdir) or whatever is appropriate
for your OS. Making directories is a function of the OS, not of R. If
you need to make a truly cross-platform solution, you might need
to check within your code what OS is being used, and call the
appropriate system statement. (I think you can do this, but have
never needed to.) That would be particularly important if you need to
specify paths.

Sarah

On 3/27/06, pau carre [EMAIL PROTECTED] wrote:

 Hello, I am trying to create directories with R. I would like R to
 create directories because it is platform independent. I tried using
 file() and searching in R Data Import/Export but I did not succeed.
 I think it must be some function since exists the unlink to remove
 directories (and files).

 Pau

 __
 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




--
Sarah Goslee
http://www.stringpage.com

[[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] seq(2,5,-2) not an error but NULL

2006-03-27 Thread Duncan Murdoch
On 3/27/2006 4:41 AM, Christian Hoffmann wrote:
 Hi,
 
 This may belong more to r-develop, but general discussion may be useful 
 (for the how many-th time ?)
 
 seq(2,5,-2)
 seq(5,2,2)
 
 both result in
 
 Error in seq.default(2, 5, -2) : wrong sign in 'by' argument
 
 But often, if not always, mathematicians and programmers want a 
 behaviour e.g. in for loops, where this statement results in an empty 
 statement, that is
 
 for (ii in seq(2,5,-2)) print(ii)
 
 were equivalent to
 
 for (ii in NULL) print(ii).
 
 The relevant part in seq.default is now
 
  if (n  0)
  stop(wrong sign in 'by' argument)
 
 but could be changed by option to
 
return(NULL)
 
 I think there should be an option to seq requiring this behaviour, or a 
 specific function, may be even a special operator, e.g. %;%:
 
 3;5 resulting in NULL.
 
 What do you think?

If you want optional behaviour, the easiest way is to write your own 
wrapper function.  E.g.

emptyseq - function(from, to, by) {
   if ((to-from)*by  0) return(NULL)
   else return(seq(from, to, by))
}

I don't think this is a desirable default, though.  We already have a 
special way to handle the most common case, i.e.

seq(1, length(x), 1)

should be written as

seq(along=x))

to handle the length(x) == 0 case the way you're requesting.

But I'm not so sure that seq(2,5,-2) should really be NULL; it looks 
much more like an error to me.  You say mathematicians and programmers 
want this behaviour, but I really can't think of any examples other than 
the one above.

As a general principle, I think it's better to throw an error on 
ambiguous or apparently erroneous code rather than always returning an 
answer.  If the code can be made unambiguous, it should be.  (R doesn't 
always follow this principle; for example, recycling of vectors of 
lengths bigger than 1 is probably an error at least as often as it's 
intended.)

Duncan Murdoch

__
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] seq(2,5,-2) not an error but NULL

2006-03-27 Thread Liaw, Andy
You should be able to do it yourself; e.g.,

my.seq - function(...) if((to - from) * by  0) NULL else seq(...)

and use that instead when you want that behavior.

Andy 

From: Christian Hoffmann
 
 Hi,
 
 This may belong more to r-develop, but general discussion may 
 be useful 
 (for the how many-th time ?)
 
 seq(2,5,-2)
 seq(5,2,2)
 
 both result in
 
 Error in seq.default(2, 5, -2) : wrong sign in 'by' argument
 
 But often, if not always, mathematicians and programmers want a 
 behaviour e.g. in for loops, where this statement results in an empty 
 statement, that is
 
 for (ii in seq(2,5,-2)) print(ii)
 
 were equivalent to
 
 for (ii in NULL) print(ii).
 
 The relevant part in seq.default is now
 
  if (n  0)
  stop(wrong sign in 'by' argument)
 
 but could be changed by option to
 
return(NULL)
 
 I think there should be an option to seq requiring this 
 behaviour, or a 
 specific function, may be even a special operator, e.g. %;%:
 
 3;5 resulting in NULL.
 
 What do you think?
 
 Christian
 -- 
 Dr. Christian W. Hoffmann,
 Swiss Federal Research Institute WSL
 Mathematics + Statistical Computing
 Zuercherstrasse 111
 CH-8903 Birmensdorf, Switzerland
 
 Tel +41-44-7392-277  (office)   -111(exchange)
 Fax +41-44-7392-215  (fax)
 [EMAIL PROTECTED] http://www.wsl.ch/staff/christian.hoffmann
 
 International Conference 5.-7.6.2006 Ekaterinburg Russia 
 Climate changes and their impact on boreal and temperate 
 forests http://ecoinf.uran.ru/conference/
 
 __
 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] seq(2,5,-2) not an error but NULL

2006-03-27 Thread Prof Brian Ripley
On Mon, 27 Mar 2006, Christian Hoffmann wrote:

 Hi,

 This may belong more to r-develop, but general discussion may be useful
 (for the how many-th time ?)

The place for general discussion of changes to R is the R-devel list.
There is almost no scope to change things like this, as there is so much 
existing code which relies on it (and it is also compatible with S).

 seq(2,5,-2)
 seq(5,2,2)

 both result in

 Error in seq.default(2, 5, -2) : wrong sign in 'by' argument

 But often, if not always, mathematicians and programmers want a
 behaviour e.g. in for loops, where this statement results in an empty
 statement, that is

 for (ii in seq(2,5,-2)) print(ii)

 were equivalent to

 for (ii in NULL) print(ii).

 The relevant part in seq.default is now

 if (n  0)
 stop(wrong sign in 'by' argument)

 but could be changed by option to

   return(NULL)

Why is NULL plausible?  I would think integer(0) is more likely, but if 
you think this should not be an error, then another plausible 
interpretation is the intersection of {2 ... 5} and {2, 0, -2, ...}, that 
is {2}.  (The only language I can think of with a close analogue is 
Fortran DO loops, where DO 10 I=2,5,-2 used to be {2} and was defined in 
F77 to be empty, so I don't think the interpretation is unambiguous.)

[As an aside, one could argue that 'for (ii in NULL)' should be an error, 
since the help page says

  seq: An expression evaluating to a vector (including a list).

and NULL is not a vector.]

 I think there should be an option to seq requiring this behaviour, or a
 specific function, may be even a special operator, e.g. %;%:

 3;5 resulting in NULL.

 What do you think?

It cannot be the default, and if you need it, why not write your own 
function to do it?  S has survived ca 18 years without it, so the need 
cannot be overwhelming.

-- 
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] How to create a directoy with R

2006-03-27 Thread Prof Brian Ripley
On Mon, 27 Mar 2006, pau carre wrote:

 Hello, I am trying to create directories with R. I would like R to
 create directories because it is platform independent. I tried using
 file() and searching in R Data Import/Export but I did not succeed.
 I think it must be some function since exists the unlink to remove
 directories (and files).

?dir.create

[help.search(directory) leads to

files(base) File and Directory Manipulation

on which page it is documented.]

-- 
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] How to create a directoy with R

2006-03-27 Thread Liaw, Andy
help.search(directory) would have given you:


R.home(base)Return the R Home Directory
files(base) File and Directory Manipulation
getwd(base) Get or Set Working Directory
list.files(base)List the Files in a Directory/Folder
unlink(base)Delete Files and Directories

and then ?files would have told you there's dir.create() and how to use it.

Andy


From: pau carre
 
 Hello, I am trying to create directories with R. I would like 
 R to create directories because it is platform independent. I 
 tried using
 file() and searching in R Data Import/Export but I did not 
 succeed. I think it must be some function since exists the 
 unlink to remove directories (and files).
 
 Pau
 
 __
 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] glmmML

2006-03-27 Thread Szentirmai Istvan
Dear R Users,

I'm looking for a similar function as step() or drop1() for glmmML models, 
but couldn't yet find any. I would appreciate if anyone could help me find 
such a function.

Thanks,
Istvan

__
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] apply(ing) to sum subset of a vector

2006-03-27 Thread Jacques VESLOT
apply(cbind(from,to), 1, function(x) sum(g[x[1]:x[2]]))

Fred J. a écrit :

Dear R users
 
 I am trying to sum selective elements of a vector but my solution
 is not cutting it.
 
 Example:
  g - 1:5;
 
  from - 1:3;
  to - 3:5;
 from to
 1   3
 2   4
 3   5
 
 so I expect 3 sums from g
 1+2+3  that is 1 to 3 of g
 2+3+4  that is 2 to 4 of g
 3+4+5  that is 3 to 5 of g
 
 my solution will not work.
 sum.em - function(g, c1, c2) sum(g[c1:c2])
 apply(g, 1, sum.em, ...) I don't think so because apply is not
 aware of the from and to. and if I f - list(g, from, to) that will not fit 
 with the second arg of apply.
 
 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

  


__
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] Plotting with date

2006-03-27 Thread pierre clauss
thanks a lot ! it runs well now
  Pierre Clauss

Gabor Grothendieck [EMAIL PROTECTED] a écrit :
  Make sure your x axis variable really is of class Date: class(x)

plot(Sys.Date() + 0:99, 1:100)

See ?str ?class ?as.Date, ?axis.Date and the help desk article in
R News 4/1 on dates for more info.

On 3/24/06, pierre clauss 
wrote:
 Hello !
 I need your help for plotting with date in the x-axis.
 I do not manage to plot temporal curves with the label of Date (01/01/1990, 
 etc) in the x-axis.
 Thanks a lot for your answer !
 Pierre Clauss.



 -

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

Re: [R] A plotting question - how to get error bars?

2006-03-27 Thread Sean Davis



On 3/27/06 6:55 AM, [EMAIL PROTECTED] [EMAIL PROTECTED] wrote:

 Dear R list,
 
 Can anyone help with a plotting question? I'm trying to display some data
 on a plot and I've almost got the format I need (see code below), but 2
 things I can't get:
 
 1. How to get Jan,Feb,Mar on the x=axis instead of 1:3?

First, do your plot with (..., axes=F).  Then, look at the help for axis()
to put the axes on the plot.

 2. How to get Ts on the end of my error bars like you have in standard
 scientific plots?

RSiteSearch('error bars') several answers that might be of interest.

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] step() for glmmML

2006-03-27 Thread Szentirmai Istvan
Dear R Users,

I'm looking for a similar function as step() or drop1() for glmmML models,
but couldn't yet find any. I would appreciate if anyone could help me find
such a function.

Thanks,
Istvan


- Original Message - 
From: [EMAIL PROTECTED]
To: r-help@stat.math.ethz.ch
Sent: Monday, March 27, 2006 1:55 PM
Subject: [R] A plotting question - how to get error bars?


 Dear R list,

 Can anyone help with a plotting question? I'm trying to display some data
 on a plot and I've almost got the format I need (see code below), but 2
 things I can't get:

 1. How to get Jan,Feb,Mar on the x=axis instead of 1:3?
 2. How to get Ts on the end of my error bars like you have in standard
 scientific plots?

 Any comments gratefully received!

 Thanks,
 Toby

 xvals=1:3 #couldn't get it to be Jan, Feb, Mar on the x-axis
 rgrT1=c(10,20,30)
 errbarsT1lo=c(0.5830952,0.3741657,0.8944272)
 errbarsT1up=errbarsT1lo
 rgrT2=c(25,30,35)
 errbarsT2lo=c(1.356466,3.535534,1.140175)
 errbarsT2up=errbarsT2lo
 minx=min(xvals);maxx=max(xvals)
 miny=min(rgrT1-errbarsT1lo,rgrT2-errbarsT2lo);maxy=max(rgrT1+errbarsT1up,rgrT2+errbarsT2up)
 plot(x=0,y=0,type=n,xlim=c(minx,maxx),ylim=c(miny,maxy),lab=c(2,20,0),bty=l,xlab=month,ylab=Relative
 Growth Rate)
 points(x=xvals,y=rgrT1,pch=21)
 symbols(x=xvals,y=rgrT1,boxplots=cbind(0,0,errbarsT1lo,errbarsT1up,0.5),inches=FALSE,add=TRUE)
 #symbols does the error bars, but without the Ts at the end. The
 boxplot command does the Ts, but you can't have them without the box in
 the middle (and you can't have different symbols for points either)
 lines(x=xvals,y=rgrT1,lty=21)
 points(x=xvals,y=rgrT2,pch=24)
 symbols(x=xvals,y=rgrT2,boxplots=cbind(0,0,errbarsT2lo,errbarsT2up,0.5),inches=FALSE,add=TRUE)
 lines(x=xvals,y=rgrT2,lty=24)
 legend(x=right,c(Treatment 1,Treatment 2),pch=c(21,24))

 __
 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 create a directoy with R

2006-03-27 Thread Philippe Grosjean
See ?dir.create, and take care at the 'recursive' argument in case you 
have to create several subdir levels at once.
Best,

Philippe Grosjean

Sarah Goslee wrote:
 I think you need to use system(mkdir) or whatever is appropriate
 for your OS. Making directories is a function of the OS, not of R. If
 you need to make a truly cross-platform solution, you might need
 to check within your code what OS is being used, and call the
 appropriate system statement. (I think you can do this, but have
 never needed to.) That would be particularly important if you need to
 specify paths.
 
 Sarah
 
 On 3/27/06, pau carre [EMAIL PROTECTED] wrote:
 
Hello, I am trying to create directories with R. I would like R to
create directories because it is platform independent. I tried using
file() and searching in R Data Import/Export but I did not succeed.
I think it must be some function since exists the unlink to remove
directories (and files).

Pau

__
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

 
 
 
 
 --
 Sarah Goslee
 http://www.stringpage.com
 
   [[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] seq(2,5,-2) not an error but NULL

2006-03-27 Thread Robin Hankin
Hi.

seq() is a complex beast indeed.   'by' being the wrong
sign is a special case of the behaviour seen in the
following code snippets, the first of which is correctly
rejected by seq(), the second of which should arguably
return a three element complex vector.


  seq(from=1,to=3,by=1+1i)
Error in n  0 : invalid comparison with complex values

  seq(from=1,to=4+3i,by=1+1i)
Error in n  0 : invalid comparison with complex values

best wishes

Robin



On 27 Mar 2006, at 13:23, Duncan Murdoch wrote:

 On 3/27/2006 4:41 AM, Christian Hoffmann wrote:
 Hi,

 This may belong more to r-develop, but general discussion may be  
 useful
 (for the how many-th time ?)

 seq(2,5,-2)
 seq(5,2,2)

 both result in

 Error in seq.default(2, 5, -2) : wrong sign in 'by' argument

 But often, if not always, mathematicians and programmers want a
 behaviour e.g. in for loops, where this statement results in an empty
 statement, that is

 for (ii in seq(2,5,-2)) print(ii)

 were equivalent to

 for (ii in NULL) print(ii).

 The relevant part in seq.default is now

  if (n  0)
  stop(wrong sign in 'by' argument)

 but could be changed by option to

return(NULL)

 I think there should be an option to seq requiring this behaviour,  
 or a
 specific function, may be even a special operator, e.g. %;%:

 3;5 resulting in NULL.

 What do you think?

 If you want optional behaviour, the easiest way is to write your own
 wrapper function.  E.g.

 emptyseq - function(from, to, by) {
if ((to-from)*by  0) return(NULL)
else return(seq(from, to, by))
 }

 I don't think this is a desirable default, though.  We already have a
 special way to handle the most common case, i.e.

 seq(1, length(x), 1)

 should be written as

 seq(along=x))

 to handle the length(x) == 0 case the way you're requesting.

 But I'm not so sure that seq(2,5,-2) should really be NULL; it looks
 much more like an error to me.  You say mathematicians and programmers
 want this behaviour, but I really can't think of any examples other  
 than
 the one above.

 As a general principle, I think it's better to throw an error on
 ambiguous or apparently erroneous code rather than always returning an
 answer.  If the code can be made unambiguous, it should be.  (R  
 doesn't
 always follow this principle; for example, recycling of vectors of
 lengths bigger than 1 is probably an error at least as often as it's
 intended.)

 Duncan Murdoch

 __
 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

--
Robin Hankin
Uncertainty Analyst
National Oceanography Centre, Southampton
European Way, Southampton SO14 3ZH, UK
  tel  023-8059-7743

__
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] seq(2,5,-2) not an error but NULL

2006-03-27 Thread Duncan Murdoch
On 3/27/2006 8:28 AM, Robin Hankin wrote:
 Hi.
 
 seq() is a complex beast indeed.   'by' being the wrong
 sign is a special case of the behaviour seen in the
 following code snippets, the first of which is correctly
 rejected by seq(), the second of which should arguably
 return a three element complex vector.
 
 
   seq(from=1,to=3,by=1+1i)
 Error in n  0 : invalid comparison with complex values
 
   seq(from=1,to=4+3i,by=1+1i)
 Error in n  0 : invalid comparison with complex values

I don't think seq() could reasonably be expected to handle to and by 
arguments with complex values.  Trying to divide the (to-from) 
difference by (by) to find how many steps to take would usually result 
in enough rounding error that the result wouldn't be real-valued.  It's 
enough of a miracle that it correctly handles

seq(from=1, by=1+1i, len=4)

Duncan Murdoch


 
 best wishes
 
 Robin
 
 
 
 On 27 Mar 2006, at 13:23, Duncan Murdoch wrote:
 
 On 3/27/2006 4:41 AM, Christian Hoffmann wrote:
 Hi,

 This may belong more to r-develop, but general discussion may be  
 useful
 (for the how many-th time ?)

 seq(2,5,-2)
 seq(5,2,2)

 both result in

 Error in seq.default(2, 5, -2) : wrong sign in 'by' argument

 But often, if not always, mathematicians and programmers want a
 behaviour e.g. in for loops, where this statement results in an empty
 statement, that is

 for (ii in seq(2,5,-2)) print(ii)

 were equivalent to

 for (ii in NULL) print(ii).

 The relevant part in seq.default is now

  if (n  0)
  stop(wrong sign in 'by' argument)

 but could be changed by option to

return(NULL)

 I think there should be an option to seq requiring this behaviour,  
 or a
 specific function, may be even a special operator, e.g. %;%:

 3;5 resulting in NULL.

 What do you think?

 If you want optional behaviour, the easiest way is to write your own
 wrapper function.  E.g.

 emptyseq - function(from, to, by) {
if ((to-from)*by  0) return(NULL)
else return(seq(from, to, by))
 }

 I don't think this is a desirable default, though.  We already have a
 special way to handle the most common case, i.e.

 seq(1, length(x), 1)

 should be written as

 seq(along=x))

 to handle the length(x) == 0 case the way you're requesting.

 But I'm not so sure that seq(2,5,-2) should really be NULL; it looks
 much more like an error to me.  You say mathematicians and programmers
 want this behaviour, but I really can't think of any examples other  
 than
 the one above.

 As a general principle, I think it's better to throw an error on
 ambiguous or apparently erroneous code rather than always returning an
 answer.  If the code can be made unambiguous, it should be.  (R  
 doesn't
 always follow this principle; for example, recycling of vectors of
 lengths bigger than 1 is probably an error at least as often as it's
 intended.)

 Duncan Murdoch

 __
 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
 
 --
 Robin Hankin
 Uncertainty Analyst
 National Oceanography Centre, Southampton
 European Way, Southampton SO14 3ZH, UK
   tel  023-8059-7743
 
 


__
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] seq(2,5,-2) not an error but NULL

2006-03-27 Thread Robin Hankin
Hi Duncan et al



 I don't think seq() could reasonably be expected to handle to and  
 by arguments with complex values.  Trying to divide the (to-from)  
 difference by (by) to find how many steps to take would usually  
 result in enough rounding error that the result wouldn't be real- 
 valued.  It's enough of a miracle that it correctly handles

 seq(from=1, by=1+1i, len=4)

 Duncan Murdoch



well it depends on your definition of miracles, but I wouldn't say

1 + (0:3) * (1+1i)

is particularly miraculous ;-\


best wishes


rksh

--
Robin Hankin
Uncertainty Analyst
National Oceanography Centre, Southampton
European Way, Southampton SO14 3ZH, UK
  tel  023-8059-7743

__
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] Date in dataframe manipulation

2006-03-27 Thread Dan Chan
Thank you Marc and Don's help, especially Marc's.  

Output- subset(FireDataAppling, select = c(STARTDATE, County, TOTAL,
CAUSE))
Worked! 
STARTDATE IS a factor and I used the following command to get the
-mm-dd format of the date
Output$Date- as.POSIXct(Output$STARTDATE)

Thank you! 

Daniel Chan

-Original Message-
From: Marc Schwartz (via MN) [mailto:[EMAIL PROTECTED] 
Sent: Friday, March 24, 2006 9:22 PM
To: Dan Chan
Cc: r-help@stat.math.ethz.ch
Subject: Re: [R] Date in dataframe manipulation

On Fri, 2006-03-24 at 15:29 -0500, Dan Chan wrote:
 Hi,
 
 I have a dataframe with many columns, including date and I want to
keep
 only a few of the columns including date column.
 
 I used the following command: 
 with(FireDataAppling, cbind(STARTDATE, County, TOTAL, CAUSE)
 
 It works, but the date becomes days from Jan 1, 2001.  
 
 FireDataAppling$STARTDATE[1] gives
 [1] 2001-01-04 00:00:00  
 1703 Levels: .

This output suggests that STARTDATE is a factor, rather than a Date
related data type. Did you read this data in via one of the read.table()
family of functions? If these values are quoted character fields in the
imported text file, they will be converted to factors by default.

 After the cbind command, the entry becomes a 4.  
 
 I want to get 2001-01-04.  What command should I use?  
 
 Thank you. 

You might want to review the Note section in ?cbind, relative to the
result of cbind()ing vectors of differing data types. By using with(),
you are effectively taking the data frame columns as individual vectors
and the resultant _matrix_ will be coerced to a single data type, in
this case, presumably numeric. I am guessing that 'County' and 'CAUSE'
are also factors, whereas 'TOTAL' is numeric.

Using str(FireDataAppling) will give you some insight into the structure
of your data frame.

The '4' that you are getting is the factor level numeric code for the
entry above, not the number of days since Jan 1, 2001, which is not a
default 'origin' date in R. Jan 1, 1970 is.

You might want to look at ?factor for more insight here.

If you want to retain only a _subset_ of the columns in a data frame,
use the subset() function:

  subset(FireDataAppling, select = c(STARTDATE, County, TOTAL, CAUSE))

This will return a data frame and retain the original data types. If you
want to then perform actual Date based operations on those values, take
a look at ?DateTimeClasses, paying attention to the See Also section
relative to associated functions.

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


Re: [R] glmmML

2006-03-27 Thread Prof Brian Ripley
yOn Mon, 27 Mar 2006, Szentirmai Istvan wrote:

 Dear R Users,

 I'm looking for a similar function as step() or drop1() for glmmML models,
 but couldn't yet find any. I would appreciate if anyone could help me find
 such a function.

If you read the help for those functions, you will see that 'all' you need 
is to specify an extractAIC() function for such models.

However, I don't think even defining AIC for them is that easy (what do 
you do about zero variance parameters, for example?).

-- 
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] Differential Equations

2006-03-27 Thread Dominik Heinzmann
Dear R-community

My ODE problems looks as follows:
(1) dA/dt = u*A - v*B
(2) dB/dt = v*B - u*A

where u is a constant, and v=k*t (k=constant, t=time)

Does anybody knows a good function/procedure of solving? Should one 
involve the equation (3) dv/dt = k?
Thanks for your support.


-- 
Dominik Heinzmann
Master of Science in Mathematics, EPFL
Ph.D. student in Biostatistics
Institute of Mathematics
University of Zurich

__
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] A plotting question - how to get error bars?

2006-03-27 Thread Gabor Grothendieck
See plotCI in package gplots.

For dates you can make use of the builtin vector month.abb

plot(1:3, 11:13, xaxt = n)
axis(1, 1:3, month.abb[1:3])

or use Date class:

xvals - seq(as.Date(2006-01-01), length = 3, by = month)
plot(xvals, 1:3)

or with specific control over x axis:

xvals - seq(as.Date(2006-01-01), length = 3, by = month)
plot(xvals, 1:3, xaxt = n)
axis(1, xvals, format(xvals, %b))

More on dates is in the Help Desk article on R News 4/1.


As an aside, note that xlim=range(xvals) is a bit more compact..


On 3/27/06, Sean Davis [EMAIL PROTECTED] wrote:



 On 3/27/06 6:55 AM, [EMAIL PROTECTED] [EMAIL PROTECTED] wrote:

  Dear R list,
 
  Can anyone help with a plotting question? I'm trying to display some data
  on a plot and I've almost got the format I need (see code below), but 2
  things I can't get:
 
  1. How to get Jan,Feb,Mar on the x=axis instead of 1:3?

 First, do your plot with (..., axes=F).  Then, look at the help for axis()
 to put the axes on the plot.

  2. How to get Ts on the end of my error bars like you have in standard
  scientific plots?

 RSiteSearch('error bars') several answers that might be of interest.

 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


[R] Graded Response Model Simulation (SAS code conversion)

2006-03-27 Thread Me
I have used R a lot in the past, but never for simulation.  I have a code in 
SAS for the Graded Response Model (GRM), also known as Samejima's model.  This 
code simulates an ordinal response, provided item characteristics (A=item 
discrimination, BB(G) are thresholds between various categorical responses).  
It is a macro file.  I am thinking that I can write this as a function, and 
call it up inside a simulation code.  Here is the SAS code:
   
  %MACRO GRGEN;
  DO G=1 TO NCAT-1;
   Z=EXP(A*(THETA-BB(G))); PS(G)=Z/(1+Z);
  END;
  PP(1)=1-PS(1); PP(NCAT)=PS(NCAT-1);
  DO G=2 TO NCAT-1;
PP(G)=PS(G-1)-PS(G);
  END;
X=RANUNI(-1);
SUMP=0; R(J)=1;
DO K=1 TO NCAT-1;
 SUMP=SUMP+PP(K);
IF XSUMP THEN R(J)=K+1;
END;
  %MEND GRGEN;
   
  Now, I am totally unfamiliar to simulation in R.  So does anyone have a good 
reference I could go to convert this?  Or have any suggestions for how to 
convert it to R?
   
  My biggest problem is all the loops inside this program.  In particular, how 
to setup the updating of R(J).
   
  It seems if I built a function for this, I need the item parameters A and 
BB's (possibly the NCAT).
   
  Any suggestions?
   
  Thanks for any help or info.
  Keith Yang
  University of Tennessee


-

[[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] comparing AIC values of models with transformed, untransformed, and weighted variables

2006-03-27 Thread Prof Brian Ripley
Two comments:

1) The log-likelihood and hence AIC for a model for log X are not 
comparable with those of a model for X.  You need to make an additive 
adjustment when you transform: it is quite easy to work out what from the 
definitions.

2) The AIC given by glm() for weighted models was wrong in R  2.3.0 
alpha.  I am not sure why you are using a glm for what appears to be a 
least-squares fit: use lm() instead (or try 2.3.0 alpha).

On Wed, 15 Mar 2006, Patrick Baker wrote:

 Hi there, I have a question regarding model comparisons that seems simple 
 enough but to which I cannot find an answer. I am interested in developing a 
 predictive model relating some measure of a tree's stem to the total leaf 
 area (TLA) of the tree. Predictor variables might include, for example, the 
 total cross-sectional area of the tree (commonly referred to as basal area) 
 or the amount of sapwood area (SA) (which represents the amount of wood 
 involved in active transport of water up the tree to the leaves). A variety 
 of people have developed these models for a variety of tree species in a 
 variety of places around the world. Perhaps not surprisingly, different 
 studies have used different model forms in analyzing their data. I am 
 interested in comparing the range of models that have been previously used 
 (some of which are theoretically derived, others of which are empirically 
 driven) using a data set that I have collected (for yet another species in 
 yet another place). To compare the different model forms I had intended to 
 use the AIC. However, I have found, again perhaps not surprisingly, that when 
 I use log-transformed data, the AIC is substantially lower for a given 
 predictor variable. If I use a weighted glm the same issue arises. For 
 example, using BA vs TLA the (rounded) AIC values are  275 for a linear 
 model, 30 for a log-log model, and 8 for a glm weighted by 1/BA. I don't 
 believe that these vast differences reflect a major improvement in the model, 
 but rather the scaling of the variables by transformation or weighting. What 
 I'd like to get some advice or insight on is whether there is an appropriate 
 way to rescale the AIC values to permit  comparisons across these models. Any 
 suggestions would be very welcome. Cheers, Patrick Baker


-- 
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] Differential Equations

2006-03-27 Thread Peter Dalgaard
Dominik Heinzmann [EMAIL PROTECTED] writes:

 Dear R-community
 
 My ODE problems looks as follows:
 (1) dA/dt = u*A - v*B
 (2) dB/dt = v*B - u*A
 
 where u is a constant, and v=k*t (k=constant, t=time)
 
 Does anybody knows a good function/procedure of solving? Should one 
 involve the equation (3) dv/dt = k?
 Thanks for your support.

You probably need to look into the odesolve package. The v*B terms
make it nonlinear, so I wouldn't expect an analytic solution.

-- 
   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] useR! 2006 program online

2006-03-27 Thread Achim Zeileis
Dear useRs,

we are happy to inform you that the program for the 2nd R user conference
useR! 2006 is now available online from the conference Web page at
  http://www.R-project.org/useR-2006/program.html

We would like to thank the useR community for submitting so many
interesting abstracts about an astonishing variety of R applications. We
are looking forward to an exciting and diversified conference!

The useR! organization team and program committee

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

__
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] X11, fonts, R-2.0.1, R-2.2.1 and R-devel

2006-03-27 Thread Prof Brian Ripley
The issue here is that you are using a UTF-8 locale (you sent this message 
in UTF-8), and you need appropriately encoded X11 fonts.  R-2.0.1 did not 
support UTF-8, and so you got incorrect output for non-ASCII characters.


It *is* an X11 installation/fontpath problem.

On Thu, 23 Mar 2006, Xavier Fernández i Marín wrote:


Hello,

I am having some problems with the X11 display in a gentoo linux laptop
with R compiled manually.
(https://stat.ethz.ch/pipermail/r-help/2006-March/089701.html)


Whether I can open the X11 device and use it when I am using 'ion' as a
window manager, I can't open it using 'gnome', due to a problem related to
fonts:
-8---
Error in X11() : could not find any X11 fonts
Check that the Font Path is correct.
-8---

I have tried and compiled R-2.0.1 and _it works_.
With the latest stable version of R it does not work.
And with the latest development (22 march 06) it does not work, neither.


It is not due to the Xorg installation, because the display is opened when
using other window manager different from gnome (and other version of R)

It is not something related to the compilation options, for the same
reason.

So my last option is that it seems to be a problem with R, X11 and gnome,
specifically.

Any ideas or suggestions?


I have googled and somebody in the FreeBSD lists talked about more or less
the same problems, but it seems that without success:
http://www.archivum.info/[EMAIL PROTECTED]/2005-06/msg00313.html
http://www.archivum.info/[EMAIL PROTECTED]/2005-06/msg00307.html

I paste the sessions and the errors, if it helps:
-8---
R-2.0.1 $ ./bin/R

R : Copyright 2004, The R Foundation for Statistical Computing
Version 2.0.1  (2004-11-15), ISBN 3-900051-07-0

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for a HTML browser interface to help.
Type 'q()' to quit R.

WARNING: UTF-8 locales are not currently supported


X11()
q()

Save workspace image? [y/n/c]: n
-8---

Works. The display is opened.

-8---
R-2.2.1 $ ./bin/R

R : Copyright 2005, The R Foundation for Statistical Computing
Version 2.2.1  (2005-12-20 r36812) ISBN 3-900051-07-0

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.


X11()

Error in X11() : could not find any X11 fonts
Check that the Font Path is correct.

q()

-8---

Doesn't work.


-8---
deu R-devel $ ./bin/R

R : Copyright 2006, The R Foundation for Statistical Computing
Version 2.3.0 Under development (unstable) (2006-03-22 r37566)
ISBN 3-900051-07-0

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.


X11()

Error in X11(display, width, height, pointsize, if (is.null(gamma)) 1 else
gamma,  :
   invalid 'width' or 'height'

x - rnorm(50)
y - rnorm(50)
plot(x,y)

Error in X11(display, width, height, pointsize, if (is.null(gamma)) 1 else
gamma,  :
   invalid 'width' or 'height'

X11(width=200, height=200)

Error in X11(width = 200, height = 200) : could not find any X11 fonts
Check that the Font Path is correct.

-8---





--
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] Arrays of functions or results of functions.

2006-03-27 Thread Gary Mallard
The general problem I am trying to solve is to determine if a series of 
subsets of data can be described with a single regression slope.  This 
involves fitting the data to each subset, calculating a joint slope 
followed by F tests to show that the variances are equal the final slope is 
valid.
The data for is characterized by a parameter PC (for the 4 - in this case) 
sets and a dependent variable RI and an independent variable ROH.
The data are contained in a variable joint.
The joint has been attached and has RI, ROH and PC for each element.
The following gives the initial results:
Mline-lm(RI[PC==1]~ROH[PC==1])
Eline-lm(RI[PC==2]~ROH[PC==2])
Iline-lm(RI[PC==3]~ROH[PC==3])
Pline-lm(RI[PC==4]~ROH[PC==4])


joint_reduced - joint;
for(i in 1:4) {

joint_reduced$RI[joint_reduced$PC==i]-joint$RI[joint$PC==i]-mean(joint$RI[joint$PC==1]}
AllLine-lm(joint_reduced$RI~joint_reduced$ROH);

Now the statistics from AllLine can be compared with each of the individual 
statistics.

NOW THE QUESTION:
 From a lot of point of view it would be useful to have a parameter 
generated by
for (i in 1:4){ Xline[i]=lm(RI[PC==i]~ROH[PC==i])}
And now all of the work of comparison can be done with calls to Xline[i] 
rather than having to  work individually with Mline, Eline, etc.

This appears to be impossible.  The constructor for Xline[i] is not 
automatic (as it is for Mline, etc) noted above.  I cannot determine how to 
construct the Xline[i] object so that this kind of process can be 
generalized.  Is it possible?  Is there another way to set us such tests of 
multiple line linearity that is already in a package?

Comments or pointers would be appreciated.
Gary

__
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] Glm poisson

2006-03-27 Thread Guenther, Cameron
Hello,
I am using the glm model with a poisson distribution.  The model runs
just fine but when I try to get the null deviance for the model of the
null degrees of freedom I get the following errors:

 null.deviance(pAmeir_1)
Error: couldn't find function null.deviance

 df.null(pAmeir_1)
Error: couldn't find function df.null

When I do:

 names(pAmeir_1)
 [1] coefficients  residuals fitted.values
 [4] effects   R rank 
 [7] qrfamilylinear.predictors
[10] deviance  aic   null.deviance
[13] iter  weights   prior.weights
[16] df.residual   df.null   y
[19] converged boundary  model
[22] call  formula   terms
[25] data  offsetcontrol  
[28] methodcontrasts xlevels   

It tells me that the values are there but I can not access them.

Any comments would be greatly appreciated.

Thank you

Cameron Guenther, Ph.D. 
Associate Research Scientist
FWC/FWRI, Marine Fisheries Research
100 8th Avenue S.E.
St. Petersburg, FL 33701
(727)896-8626 Ext. 4305
[EMAIL PROTECTED]

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


Re: [R] Differential Equations

2006-03-27 Thread Ravi Varadhan
Dominik,

Adding (1) and (2) yields,

A(t) + B(t) = constant = A(0) + B(0) = c

So, plug in B = c - A in (1) and solve for A.  This should be an easy
solution.

Hope this is helpful,
Ravi.

 -Original Message-
 From: [EMAIL PROTECTED] [mailto:r-help-
 [EMAIL PROTECTED] On Behalf Of Peter Dalgaard
 Sent: Monday, March 27, 2006 11:37 AM
 To: Dominik Heinzmann
 Cc: r-help@stat.math.ethz.ch
 Subject: Re: [R] Differential Equations
 
 Dominik Heinzmann [EMAIL PROTECTED] writes:
 
  Dear R-community
 
  My ODE problems looks as follows:
  (1) dA/dt = u*A - v*B
  (2) dB/dt = v*B - u*A
 
  where u is a constant, and v=k*t (k=constant, t=time)
 
  Does anybody knows a good function/procedure of solving? Should one
  involve the equation (3) dv/dt = k?
  Thanks for your support.
 
 You probably need to look into the odesolve package. The v*B terms
 make it nonlinear, so I wouldn't expect an analytic solution.
 
 --
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-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] Arrays of functions or results of functions.

2006-03-27 Thread Gabor Grothendieck
Is your question how to run a regression with separate slopes
and then with one slope and then to complare them?  If that
is it here is an example:

# test data - ind is factor which defines the groups
set.seed(1)
y1 - 10 + 20 * seq(100) + rnorm(100)
y2 - -200 + 35 * seq(100) + rnorm(100)
yy - pmax(y1, y2)
ind - factor(y1  y2)

# lm.N has N slopes
lm.1 - lm(yy ~ ind + seq(100) - 1)
lm.2 - lm(yy ~ ind/seq(100)-1)
lm.1
lm.2
anova(lm.1, lm.2)

On 3/27/06, Gary Mallard [EMAIL PROTECTED] wrote:
 The general problem I am trying to solve is to determine if a series of
 subsets of data can be described with a single regression slope.  This
 involves fitting the data to each subset, calculating a joint slope
 followed by F tests to show that the variances are equal the final slope is
 valid.
 The data for is characterized by a parameter PC (for the 4 - in this case)
 sets and a dependent variable RI and an independent variable ROH.
 The data are contained in a variable joint.
 The joint has been attached and has RI, ROH and PC for each element.
 The following gives the initial results:
 Mline-lm(RI[PC==1]~ROH[PC==1])
 Eline-lm(RI[PC==2]~ROH[PC==2])
 Iline-lm(RI[PC==3]~ROH[PC==3])
 Pline-lm(RI[PC==4]~ROH[PC==4])


 joint_reduced - joint;
 for(i in 1:4) {

 joint_reduced$RI[joint_reduced$PC==i]-joint$RI[joint$PC==i]-mean(joint$RI[joint$PC==1]}
 AllLine-lm(joint_reduced$RI~joint_reduced$ROH);

 Now the statistics from AllLine can be compared with each of the individual
 statistics.

 NOW THE QUESTION:
  From a lot of point of view it would be useful to have a parameter
 generated by
 for (i in 1:4){ Xline[i]=lm(RI[PC==i]~ROH[PC==i])}
 And now all of the work of comparison can be done with calls to Xline[i]
 rather than having to  work individually with Mline, Eline, etc.

 This appears to be impossible.  The constructor for Xline[i] is not
 automatic (as it is for Mline, etc) noted above.  I cannot determine how to
 construct the Xline[i] object so that this kind of process can be
 generalized.  Is it possible?  Is there another way to set us such tests of
 multiple line linearity that is already in a package?

 Comments or pointers would be appreciated.
 Gary

 __
 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] Glm poisson

2006-03-27 Thread Francisco J. Zagmutt
Why are you trying to extract the values by calling a function with the name 
of the value?  glm objects are stored as a list i.e.
str(pAmeir_1)

Hence, you can extract what you need by selecting the values on the list 
i.e.

pAmeir_1$df.null
pAmeir_1$null.deviance

Cheers

Francisco

From: Guenther, Cameron [EMAIL PROTECTED]
To: [EMAIL PROTECTED], r-help@stat.math.ethz.ch
Subject: [R] Glm poisson
Date: Mon, 27 Mar 2006 12:07:04 -0500

Hello,
I am using the glm model with a poisson distribution.  The model runs
just fine but when I try to get the null deviance for the model of the
null degrees of freedom I get the following errors:

  null.deviance(pAmeir_1)
Error: couldn't find function null.deviance

  df.null(pAmeir_1)
Error: couldn't find function df.null

When I do:

  names(pAmeir_1)
  [1] coefficients  residuals fitted.values
  [4] effects   R rank
  [7] qrfamilylinear.predictors
[10] deviance  aic   null.deviance
[13] iter  weights   prior.weights
[16] df.residual   df.null   y
[19] converged boundary  model
[22] call  formula   terms
[25] data  offsetcontrol
[28] methodcontrasts xlevels

It tells me that the values are there but I can not access them.

Any comments would be greatly appreciated.

Thank you

Cameron Guenther, Ph.D.
Associate Research Scientist
FWC/FWRI, Marine Fisheries Research
100 8th Avenue S.E.
St. Petersburg, FL 33701
(727)896-8626 Ext. 4305
[EMAIL PROTECTED]

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

__
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] Arrays of functions or results of functions.

2006-03-27 Thread Berton Gunter
If I understand you correctly, I would say that this is a standard analysis
of covariance problem that you are approaching incorrectly. You should not
be testing for equal variances, IMO. Instead,

1. Combine all your data into 3 columnsS x, y, and group= subset.

2. Model.1: y ~ x
3. Model.2: y ~ x*group.  This gives a separate slope and intercept for each
group.
4. anova(Model.1, Model.2)  to compare.

A better approach would be to use lmList() from the nlme package. Better yet
would be model the data as a mixed effect model via lme(), which assumes
that slopes and intercepts vary randomly among your subsets about their
corresponfding central values. This gets you shrinkage, which is highly
desirable. See Bates and Pinheiro MIXED EFFECTS MODELS IN S AND S-PLUS for
details.

-- Bert Gunter
Genentech Non-Clinical Statistics
South San Francisco, CA
 
The business of the statistician is to catalyze the scientific learning
process.  - George E. P. Box
 
 

 -Original Message-
 From: [EMAIL PROTECTED] 
 [mailto:[EMAIL PROTECTED] On Behalf Of Gabor 
 Grothendieck
 Sent: Monday, March 27, 2006 9:48 AM
 To: Gary Mallard
 Cc: r-help@stat.math.ethz.ch
 Subject: Re: [R] Arrays of functions or results of functions.
 
 Is your question how to run a regression with separate slopes
 and then with one slope and then to complare them?  If that
 is it here is an example:
 
 # test data - ind is factor which defines the groups
 set.seed(1)
 y1 - 10 + 20 * seq(100) + rnorm(100)
 y2 - -200 + 35 * seq(100) + rnorm(100)
 yy - pmax(y1, y2)
 ind - factor(y1  y2)
 
 # lm.N has N slopes
 lm.1 - lm(yy ~ ind + seq(100) - 1)
 lm.2 - lm(yy ~ ind/seq(100)-1)
 lm.1
 lm.2
 anova(lm.1, lm.2)
 
 On 3/27/06, Gary Mallard [EMAIL PROTECTED] wrote:
  The general problem I am trying to solve is to determine if 
 a series of
  subsets of data can be described with a single regression 
 slope.  This
  involves fitting the data to each subset, calculating a joint slope
  followed by F tests to show that the variances are equal 
 the final slope is
  valid.
  The data for is characterized by a parameter PC (for the 4 
 - in this case)
  sets and a dependent variable RI and an independent variable ROH.
  The data are contained in a variable joint.
  The joint has been attached and has RI, ROH and PC for each element.
  The following gives the initial results:
  Mline-lm(RI[PC==1]~ROH[PC==1])
  Eline-lm(RI[PC==2]~ROH[PC==2])
  Iline-lm(RI[PC==3]~ROH[PC==3])
  Pline-lm(RI[PC==4]~ROH[PC==4])
 
 
  joint_reduced - joint;
  for(i in 1:4) {
 
 joint_reduced$RI[joint_reduced$PC==i]-joint$RI[joint$PC==i]-m
 ean(joint$RI[joint$PC==1]}
  AllLine-lm(joint_reduced$RI~joint_reduced$ROH);
 
  Now the statistics from AllLine can be compared with each 
 of the individual
  statistics.
 
  NOW THE QUESTION:
   From a lot of point of view it would be useful to have a parameter
  generated by
  for (i in 1:4){ Xline[i]=lm(RI[PC==i]~ROH[PC==i])}
  And now all of the work of comparison can be done with 
 calls to Xline[i]
  rather than having to  work individually with Mline, Eline, etc.
 
  This appears to be impossible.  The constructor for Xline[i] is not
  automatic (as it is for Mline, etc) noted above.  I cannot 
 determine how to
  construct the Xline[i] object so that this kind of process can be
  generalized.  Is it possible?  Is there another way to set 
 us such tests of
  multiple line linearity that is already in a package?
 
  Comments or pointers would be appreciated.
  Gary
 
  __
  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-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] reading in multi-dimensional data from .csv

2006-03-27 Thread Werner Wernersen
Hi,

I would like to read in multi-dimensional data from a
text file, i.e. tables with more than 2 dimensions.
I have looked for a function which I can abuse for
that but haven't found anything. 

I would appreciate it a lot if somebody gave me a hint
if such functions already exist somewhere.

Thanks,
  Werner

__
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] reading in multi-dimensional data from .csv

2006-03-27 Thread jim holtman
How is the data organized?  You could 'linearize' the data (e.g., column
order) and supply the dimensions of the data that you would apply after the
data was read in.

If you had three dimensions and you separated each 2D page with a blank
line, you again could reconstruct the data.  If you have the ability to
change the program doing the output, then look at the output of 'dput' for a
multidimensional object and create a text string that looks like that and
'source' it in.


On 3/27/06, Werner Wernersen [EMAIL PROTECTED] wrote:

 Hi,

 I would like to read in multi-dimensional data from a
 text file, i.e. tables with more than 2 dimensions.
 I have looked for a function which I can abuse for
 that but haven't found anything.

 I would appreciate it a lot if somebody gave me a hint
 if such functions already exist somewhere.

 Thanks,
 Werner

 __
 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




--
Jim Holtman
Cincinnati, OH
+1 513 646 9390 (Cell)
+1 513 247 0281 (Home)

What the problem you are trying to solve?

[[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] graphing and scrolling

2006-03-27 Thread Fred J.
Dear R users
 
  graphing with plot(x) seams to work for a small length(x), when length(x)  is 
too large it seams to clutter the display, a solution would be to  display 
subsets of x at a time, yet a better way which I hope R  supports is to place a 
sliding bar on the display window to control  length(x) and thus the 
resolution, which will involve auto scaling  the y axis as well as 
automatically place a side-to-side scrolling bar  at the bottom of the chart to 
move through and display the subsets of  x  according  to the resolution. is 
there such an implementation?
 
 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


[R] DSC 2007

2006-03-27 Thread Paul Murrell
Hi

Following on from the Distributed Statistical Computing conferences in 
Vienna (1999, 2001, 2003) and the Directions in Statistical Computing 
conference in Seattle last year ...

DSC 2007, a conference on systems and environments for statistical 
computing, will take place in Auckland, New Zealand on February 15  16, 
2007.

This is just an announcement of date and location so that interested 
parties living on the far side of the planet have time to book their 
travel;  further details and a web site for the conference will be 
announced in due course.

Paul
-- 
Dr Paul Murrell
Department of Statistics
The University of Auckland
Private Bag 92019
Auckland
New Zealand
64 9 3737599 x85392
[EMAIL PROTECTED]
http://www.stat.auckland.ac.nz/~paul/

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

__
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] graphing and scrolling

2006-03-27 Thread Gregory Snow
Does the following (or some simple modification of it) do what you
want?:

library(tkrplot)

y - rnorm(1, 10, 2) + 5*sin( (1:1)/1000 )

tt - tktoplevel()
left - tclVar(1)
oldleft - tclVar(1)
right - tclVar(100)

f1 - function(){
lleft - as.numeric(tclvalue(left))
rright - as.numeric(tclvalue(right))
x - seq(lleft,rright,by=1)
plot(x,y[x], type='b',ylim=range(y))
}

img - tkrplot(tt, f1)


f2 - function(...){
ol - as.numeric(tclvalue(oldleft))
tclvalue(oldleft) - tclvalue(left)
r - as.numeric(tclvalue(right))
tclvalue(right) - as.character(r + as.numeric(...) - ol)
tkrreplot(img)
}

f3 - function(...){
tkrreplot(img)
}

f4 - function(...){
ol - as.numeric(tclvalue(oldleft))
tclvalue(left) - as.character(ol+100)
tclvalue(oldleft) - as.character(ol+100)
r - as.numeric(tclvalue(right))
tclvalue(right) - as.character(r+100)
tkrreplot(img)
}

s1 - tkscale(tt, command=f2, from=1, to=length(y),
variable=left, orient=horiz,label='left')
s2  - tkscale(tt, command=f3, from=1, to=length(y),
variable=right, orient=horiz,label='right')
b1 - tkbutton(tt, text='-', command=f4)

tkpack(img,s1,s2,b1)





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

 -Original Message-
 From: [EMAIL PROTECTED] 
 [mailto:[EMAIL PROTECTED] On Behalf Of Fred J.
 Sent: Monday, March 27, 2006 1:05 PM
 To: r-help@stat.math.ethz.ch
 Subject: [R] graphing and scrolling
 
 Dear R users
  
   graphing with plot(x) seams to work for a small length(x), 
 when length(x)  is too large it seams to clutter the display, 
 a solution would be to  display subsets of x at a time, yet a 
 better way which I hope R  supports is to place a sliding bar 
 on the display window to control  length(x) and thus the 
 resolution, which will involve auto scaling  the y axis as 
 well as automatically place a side-to-side scrolling bar  at 
 the bottom of the chart to move through and display the 
 subsets of  x  according  to the resolution. is there such an 
 implementation?
  
  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


__
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] why does lmer not give p valuses for quasibinomial family?

2006-03-27 Thread Szentirmai Istvan
Dear All,

I'm running a binomial model using the lmer() function, and I get p values 
for the parameter estimates only with family=binomial, but not with 
quasibinomial? Why is that so? I wanted to use quasibinomial family, because 
my data were overdispersed.

Thanks,
Istvan

__
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] products and polynomials in formulae

2006-03-27 Thread stephenc
Hi
 
I can do this:
 
formula = as.factor(outcome) ~ .
 
in glm and other model building functions.   I think there is a way to
get the product of the determinants (that is d1 * d2, d1 * d3, etc) and
also another way to get all the polynomials (that is like poly(d1,2)
would produce for a single determinant).
 
Can anyone tell me how you write them?
 
Stephen 

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

2006-03-27 Thread stephenc
Hi
 
Does anyone know what this means:
 
 
 glm.model = glm(formula = as.factor(nextDay) ~ ., family=binomial,
data=spi[1:1000,])
 pred - predict(glm.model, spi[1001:1250,-9], type=response)
Warning message: 
prediction from a rank-deficient fit may be misleading in:
predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type ==  
 
9 is my determinant and I still get this message even when I remove the
9.
 
Stephen
 

[[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 understanding behavior of apply vs sapply

2006-03-27 Thread Seth Falcon
Hi,

I was surprised that apply and sapply don't return the same results in
the example below.  Can someone tell me what I'm missing?


 zls - function(x) character(0)
 m - matrix(0, nrow=2, ncol=2)
 apply(m, 1, zls)
character(0)
 sapply(m, zls)
[[1]]
character(0)

[[2]]
character(0)

[[3]]
character(0)

[[4]]
character(0)






 R.version
   _  
platform   powerpc-apple-darwin8.5.0  
arch   powerpc
os darwin8.5.0
system powerpc, darwin8.5.0   
status alpha  
major  2  
minor  3.0
year   2006   
month  03 
day27 
svn rev37590  
language   R  
version.string Version 2.3.0 alpha (2006-03-27 r37590)

__
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] fixed effects

2006-03-27 Thread ivo welch
dear R wizards:

X is factor with 20,000*20=800,000 observations of 20,000 factors. 
I.e., each factor has 20 observations.  y is 800,000 normally
distributed data points.  I want to see how much R^2 the X factors can
provide.  Easy, right?

 lm ( y ~ X)
and
  aov( y ~ X)
Error: cannot allocate vector of size 3125000 Kb

is this computationally infeasible?  (I am not an expert, but
off-hand, I thought this can work as long as the X's are just factors
= fitted means.)

 help.search(fixed.effects);

fixed.effects(nlme) Extract Fixed Effects
fixed.effects.lmList(nlme)  Extract lmList Fixed Effects
lme(nlme)   Linear Mixed-Effects Models
lmeStruct(nlme) Linear Mixed-Effects Structure
nlme(nlme)  Nonlinear Mixed-Effects Models
nlmeStruct(nlme)Nonlinear Mixed-Effects Structure

Type 'help(FOO, package = PKG)' to inspect entry 'FOO(PKG) TITLE'.

ok---I want to read the fixed.effects help.  could the help system
tell me how to inspect these entries?

 help(fixed.effects)
wrong
 help.search(fixed.effects)
two entries, as above.  nothing new.
 ?fixed.effects
wrong

eventually, it dawned on me that nlme in parens was not a function
argument, but the name of the package within which fixed.effects
lives.  Suggestion:  maybe a different notation to name packages would
be more intuitive in the help system.  yes, I know it now, but other
novices may not.  even a colon instead of a () may be more intuitive
in this context.

 library(nlme); ?lme
and then
 lme(y ~ X)
Error in getGroups.data.frame(dataMix, groups) :
Invalid formula for groups


now I have to beg for some help.  ok, blatant free-riding.  the lme
docs tells me it does the Laird and Ware model, but I do not know this
model.  the only two examples given at the end of the lme help file
seem to be similar to what I just specified.  so, how do I execute a
simple fixed-effects model?  (later, I may want to add a variable Z
that is a continuous random variable.)  could someone please give me
one quick example?  help is, as always, highly appreciated.

sincerely,

/ivo welch

__
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 understanding behavior of apply vs sapply

2006-03-27 Thread Liaw, Andy
I'll give it a shot:

apply() is only a wrapper around a for loop through the requested
dimension(s).  In this case it would run zls() once for each row in m; i.e.,
twice.  The `Value' section of apply() explains what it does if the function
being applied returns vector of length 0.

sapply() is basically lapply() with nice post-processing if possible.  It
doesn't know about m being 2x2, just that it's a vector with four elements,
so it's going to call the function four times, and return a list with the
result of each of those calls.

Andy

From: Seth Falcon
 
 Hi,
 
 I was surprised that apply and sapply don't return the same 
 results in the example below.  Can someone tell me what I'm missing?
 
 
  zls - function(x) character(0)
  m - matrix(0, nrow=2, ncol=2)
  apply(m, 1, zls)
 character(0)
  sapply(m, zls)
 [[1]]
 character(0)
 
 [[2]]
 character(0)
 
 [[3]]
 character(0)
 
 [[4]]
 character(0)
 
 
 
 
 
 
  R.version
_  
 platform   powerpc-apple-darwin8.5.0  
 arch   powerpc
 os darwin8.5.0
 system powerpc, darwin8.5.0   
 status alpha  
 major  2  
 minor  3.0
 year   2006   
 month  03 
 day27 
 svn rev37590  
 language   R  
 version.string Version 2.3.0 alpha (2006-03-27 r37590)
 
 __
 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] fixed effects

2006-03-27 Thread Liaw, Andy
I guess you meant X has 20,000 levels (and 40 observations each)?  In that
case lm() will attempt to create a design matrix that's 8e5 by 2e4, and
that's unlikely to fit in the RAM.

It is very easy to compute by hand.  I'm using a smaller data size first to
check the result against summary(lm()), then the size you specified
(hopefully with more correct arithmetics...) to show that it's quite doable
even on a modest laptop:

 set.seed(1)
 y - rnorm(80)
 x - factor(rep(paste(x, 1:20, sep=), 4))
 summary(lm(y ~ x))$r.squared
[1] 0.2555885
 rsq - function(x, y) {
+ sse - sum(ave(y, x, FUN=function(x) x - mean(x))^2)
+ sstotal - var(y) * (length(y) - 1)
+ 1 - sse / sstotal
+ }
 rsq(x, y)
[1] 0.2555885
 set.seed(1)
 y - rnorm(8e5)
 x - factor(rep(paste(x, 1:2e4, sep=), 40))
 system.time(rsq(x, y))
[1] 1.99 0.03 2.06   NA   NA

Andy
 

From: ivo welch
 
 dear R wizards:
 
 X is factor with 20,000*20=800,000 observations of 20,000 factors. 
 I.e., each factor has 20 observations.  y is 800,000 normally 
 distributed data points.  I want to see how much R^2 the X 
 factors can provide.  Easy, right?
 
  lm ( y ~ X)
 and
   aov( y ~ X)
 Error: cannot allocate vector of size 3125000 Kb
 
 is this computationally infeasible?  (I am not an expert, but 
 off-hand, I thought this can work as long as the X's are just 
 factors = fitted means.)
 
  help.search(fixed.effects);
 
 fixed.effects(nlme) Extract Fixed Effects
 fixed.effects.lmList(nlme)  Extract lmList Fixed Effects
 lme(nlme)   Linear Mixed-Effects Models
 lmeStruct(nlme) Linear Mixed-Effects Structure
 nlme(nlme)  Nonlinear Mixed-Effects Models
 nlmeStruct(nlme)Nonlinear Mixed-Effects Structure
 
 Type 'help(FOO, package = PKG)' to inspect entry 'FOO(PKG) TITLE'.
 
 ok---I want to read the fixed.effects help.  could the help 
 system tell me how to inspect these entries?
 
  help(fixed.effects)
 wrong
  help.search(fixed.effects)
 two entries, as above.  nothing new.
  ?fixed.effects
 wrong
 
 eventually, it dawned on me that nlme in parens was not a 
 function argument, but the name of the package within which 
 fixed.effects lives.  Suggestion:  maybe a different notation 
 to name packages would be more intuitive in the help system.  
 yes, I know it now, but other novices may not.  even a colon 
 instead of a () may be more intuitive in this context.
 
  library(nlme); ?lme
 and then
  lme(y ~ X)
 Error in getGroups.data.frame(dataMix, groups) :
 Invalid formula for groups
 
 
 now I have to beg for some help.  ok, blatant free-riding.  
 the lme docs tells me it does the Laird and Ware model, but I 
 do not know this model.  the only two examples given at the 
 end of the lme help file seem to be similar to what I just 
 specified.  so, how do I execute a simple fixed-effects 
 model?  (later, I may want to add a variable Z that is a 
 continuous random variable.)  could someone please give me 
 one quick example?  help is, as always, highly appreciated.
 
 sincerely,
 
 /ivo welch
 
 __
 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] fixed effects

2006-03-27 Thread Liaw, Andy
Hmm... didn't read the whole post before I hit `send'...

I think you basically will have to fit the model `by hand', which is not
that hard, given the simple structure of the model(s).  The formulae for the
quantities of interests are quite straightforward and easy to code in R
(similar to the rsq function I showed).  Adding a continuous variable to the
model simply requires regressing the residuals; i.e., ave(y, x, function(x)
x - mean(x)), on the new variable z (easy for lm()), and take one more
degree of freedom from SS(total).

Andy

From: Liaw, Andy
 
 I guess you meant X has 20,000 levels (and 40 observations 
 each)?  In that case lm() will attempt to create a design 
 matrix that's 8e5 by 2e4, and that's unlikely to fit in the RAM.
 
 It is very easy to compute by hand.  I'm using a smaller data 
 size first to check the result against summary(lm()), then 
 the size you specified (hopefully with more correct 
 arithmetics...) to show that it's quite doable even on a 
 modest laptop:
 
  set.seed(1)
  y - rnorm(80)
  x - factor(rep(paste(x, 1:20, sep=), 4))
  summary(lm(y ~ x))$r.squared
 [1] 0.2555885
  rsq - function(x, y) {
 + sse - sum(ave(y, x, FUN=function(x) x - mean(x))^2)
 + sstotal - var(y) * (length(y) - 1)
 + 1 - sse / sstotal
 + }
  rsq(x, y)
 [1] 0.2555885
  set.seed(1)
  y - rnorm(8e5)
  x - factor(rep(paste(x, 1:2e4, sep=), 40)) 
 system.time(rsq(x, y))
 [1] 1.99 0.03 2.06   NA   NA
 
 Andy
  
 
 From: ivo welch
  
  dear R wizards:
  
  X is factor with 20,000*20=800,000 observations of 20,000 factors.
  I.e., each factor has 20 observations.  y is 800,000 normally 
  distributed data points.  I want to see how much R^2 the X 
  factors can provide.  Easy, right?
  
   lm ( y ~ X)
  and
aov( y ~ X)
  Error: cannot allocate vector of size 3125000 Kb
  
  is this computationally infeasible?  (I am not an expert, but
  off-hand, I thought this can work as long as the X's are just 
  factors = fitted means.)
  
   help.search(fixed.effects);
  
  fixed.effects(nlme) Extract Fixed Effects
  fixed.effects.lmList(nlme)  Extract lmList Fixed Effects
  lme(nlme)   Linear Mixed-Effects Models
  lmeStruct(nlme) Linear Mixed-Effects Structure
  nlme(nlme)  Nonlinear Mixed-Effects Models
  nlmeStruct(nlme)Nonlinear Mixed-Effects 
 Structure
  
  Type 'help(FOO, package = PKG)' to inspect entry 'FOO(PKG) TITLE'.
  
  ok---I want to read the fixed.effects help.  could the help
  system tell me how to inspect these entries?
  
   help(fixed.effects)
  wrong
   help.search(fixed.effects)
  two entries, as above.  nothing new.
   ?fixed.effects
  wrong
  
  eventually, it dawned on me that nlme in parens was not a
  function argument, but the name of the package within which 
  fixed.effects lives.  Suggestion:  maybe a different notation 
  to name packages would be more intuitive in the help system.  
  yes, I know it now, but other novices may not.  even a colon 
  instead of a () may be more intuitive in this context.
  
   library(nlme); ?lme
  and then
   lme(y ~ X)
  Error in getGroups.data.frame(dataMix, groups) :
  Invalid formula for groups
  
  
  now I have to beg for some help.  ok, blatant free-riding.
  the lme docs tells me it does the Laird and Ware model, but I 
  do not know this model.  the only two examples given at the 
  end of the lme help file seem to be similar to what I just 
  specified.  so, how do I execute a simple fixed-effects 
  model?  (later, I may want to add a variable Z that is a 
  continuous random variable.)  could someone please give me 
  one quick example?  help is, as always, highly appreciated.
  
  sincerely,
  
  /ivo welch
  
  __
  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
 
 --
 
 Notice:  This e-mail message, together with any attachments, 
 contains information of Merck  Co., Inc. (One Merck Drive, 
 Whitehouse Station, New Jersey, USA 08889), and/or its 
 affiliates (which may be known outside the United States as 
 Merck Frosst, Merck Sharp  Dohme or MSD and in Japan, as 
 Banyu) that may be confidential, proprietary copyrighted 
 and/or legally privileged. It is intended solely for the use 
 of the individual or entity named on this message.  If you 
 are not the intended recipient, and have received this 
 message in error, please notify us immediately by reply 
 e-mail and then delete it from your system.
 

Re: [R] fixed effects

2006-03-27 Thread Spencer Graves
  Have you tried the following:

  lme(y~1, random=~1|X, data=DF)

where DF = a data.frame with columns y and X.

  The authoritative reference on library(nlme) is Pinheiro and Bates 
(2000) Mixed-Effects Models in S and S-Plus (Springer).  I've learned a 
lot from Bates and from this book in particular.

  hope this helps,
  spencer graves

ivo welch wrote:

 dear R wizards:
 
 X is factor with 20,000*20=800,000 observations of 20,000 factors. 
 I.e., each factor has 20 observations.  y is 800,000 normally
 distributed data points.  I want to see how much R^2 the X factors can
 provide.  Easy, right?
 
 
lm ( y ~ X)
 
 and
 
 aov( y ~ X)
 
 Error: cannot allocate vector of size 3125000 Kb
 
 is this computationally infeasible?  (I am not an expert, but
 off-hand, I thought this can work as long as the X's are just factors
 = fitted means.)
 
 
help.search(fixed.effects);
 
 
 fixed.effects(nlme) Extract Fixed Effects
 fixed.effects.lmList(nlme)  Extract lmList Fixed Effects
 lme(nlme)   Linear Mixed-Effects Models
 lmeStruct(nlme) Linear Mixed-Effects Structure
 nlme(nlme)  Nonlinear Mixed-Effects Models
 nlmeStruct(nlme)Nonlinear Mixed-Effects Structure
 
 Type 'help(FOO, package = PKG)' to inspect entry 'FOO(PKG) TITLE'.
 
 ok---I want to read the fixed.effects help.  could the help system
 tell me how to inspect these entries?
 
 
help(fixed.effects)
 
 wrong
 
help.search(fixed.effects)
 
 two entries, as above.  nothing new.
 
?fixed.effects
 
 wrong
 
 eventually, it dawned on me that nlme in parens was not a function
 argument, but the name of the package within which fixed.effects
 lives.  Suggestion:  maybe a different notation to name packages would
 be more intuitive in the help system.  yes, I know it now, but other
 novices may not.  even a colon instead of a () may be more intuitive
 in this context.
 
 
library(nlme); ?lme
 
 and then
 
lme(y ~ X)
 
 Error in getGroups.data.frame(dataMix, groups) :
 Invalid formula for groups
 
 
 now I have to beg for some help.  ok, blatant free-riding.  the lme
 docs tells me it does the Laird and Ware model, but I do not know this
 model.  the only two examples given at the end of the lme help file
 seem to be similar to what I just specified.  so, how do I execute a
 simple fixed-effects model?  (later, I may want to add a variable Z
 that is a continuous random variable.)  could someone please give me
 one quick example?  help is, as always, highly appreciated.
 
 sincerely,
 
 /ivo welch
 
 __
 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] Having trouble with tsdiag function on a time series

2006-03-27 Thread Rafael Algara
Hello,

I'm getting the following error message when I try to run 'tsdiag' on what 
seems to be a valid time series:

 tsdiag(small)

returns:

[Error in tsdiag(small) : no applicable method for tsdiag]

where small is a little test series where I have isolated this problem (the 
original has 30-years worth of daily data)

When I print it (small), it looks like the following:

Time Series:
Start = c(1990, 1) 
End = c(1990, 31) 
Frequency = 365 
[1]  0.0 16.0 10.0  0.0  0.0  0.0  0.0  0.0  2.0  2.2  0.0  0.0  0.0  0.0  0.0  
0.0  0.0  2.5  0.1
[20]  2.4  0.0  0.0  0.0  0.0  9.6  0.0  0.0  0.0  0.0 25.2  7.0

I've found some postings in the r-bugs list which suggest that no applicable 
method can stem from a conflict across namespaces. I've tried qualifying the 
call:

 stats::tsdiag(small)

with the same result.

Many thanks in advance for any help.

Rafael
[[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] Default lag.max in ACF

2006-03-27 Thread Prof Brian Ripley
The default is taken from S-PLUS, so the reference is the S-PLUS manual.
It is pretty similar to the recommendation of Brockwell  Davis.

On Mon, 27 Mar 2006, Spencer Graves wrote:

 I don't know why the default lag.max is 10*log10(N/m) for m series. 
 The acf help page includes the following:

 Author(s):

 Original: Paul Gilbert, Martyn Plummer. Extensive modifications
 and univariate case of 'pacf' by B.D. Ripley.

 I've copied these three contributors on this email.  With luck, one 
 of them will enlighten us both.

 Best Wishes,
 spencer graves

 Rafal Stankiewicz wrote:

 Hi,
 
 The default value for lag.max in ACF implementation is 10*log10(N)
 
 There several publications recommending setting lag.max to:
 -  N/4 (Box and Jenkins, 1970; Chatfield, 1975; Anderson, 1976; Pankratz, 
 1983; Davis, 1986; etc.)
 -  sqrt(N)+10 (Cryer, 1986)
 -  20=N=40 (Brockwell and Davis)
 
 Why R uses 10*log10(N) as a default?
 
 Please, give me a  reference to a book or article where the recommendation 
 for using lag.max=10*log10(N) is proposed and explained.
 
 Thanks
 Rafal

-- 
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] Having trouble with tsdiag function on a time series

2006-03-27 Thread Prof Brian Ripley
tsdiag is for diagnostic plots from time-series fits, not time series.

On Tue, 28 Mar 2006, Rafael Algara wrote:

 Hello,

 I'm getting the following error message when I try to run 'tsdiag' on what 
 seems to be a valid time series:

 tsdiag(small)

 returns:

 [Error in tsdiag(small) : no applicable method for tsdiag]

 where small is a little test series where I have isolated this problem (the 
 original has 30-years worth of daily data)

 When I print it (small), it looks like the following:

 Time Series:
 Start = c(1990, 1)
 End = c(1990, 31)
 Frequency = 365
 [1]  0.0 16.0 10.0  0.0  0.0  0.0  0.0  0.0  2.0  2.2  0.0  0.0  0.0  0.0  
 0.0  0.0  0.0  2.5  0.1
 [20]  2.4  0.0  0.0  0.0  0.0  9.6  0.0  0.0  0.0  0.0 25.2  7.0

 I've found some postings in the r-bugs list which suggest that no 
 applicable method can stem from a conflict across namespaces. I've 
 tried qualifying the call:

Not relevant.  Try

 methods(tsdiag)
[1] tsdiag.Arima*tsdiag.arima0*   tsdiag.StructTS*

Non-visible functions are asterisked


 stats::tsdiag(small)

 with the same result.

 Many thanks in advance for any help.

 Rafael
   [[alternative HTML version deleted]]

 PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

PLEASE do, and do not send HTML mail as we ask.


-- 
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] fixed effects

2006-03-27 Thread Prof Brian Ripley
There is a very similar example worked through in section 7.2 of `S 
Programming'.

On Mon, 27 Mar 2006, Liaw, Andy wrote:

 Hmm... didn't read the whole post before I hit `send'...

 I think you basically will have to fit the model `by hand', which is not
 that hard, given the simple structure of the model(s).  The formulae for the
 quantities of interests are quite straightforward and easy to code in R
 (similar to the rsq function I showed).  Adding a continuous variable to the
 model simply requires regressing the residuals; i.e., ave(y, x, function(x)
 x - mean(x)), on the new variable z (easy for lm()), and take one more
 degree of freedom from SS(total).

 Andy

 From: Liaw, Andy

 I guess you meant X has 20,000 levels (and 40 observations
 each)?  In that case lm() will attempt to create a design
 matrix that's 8e5 by 2e4, and that's unlikely to fit in the RAM.

 It is very easy to compute by hand.  I'm using a smaller data
 size first to check the result against summary(lm()), then
 the size you specified (hopefully with more correct
 arithmetics...) to show that it's quite doable even on a
 modest laptop:

 set.seed(1)
 y - rnorm(80)
 x - factor(rep(paste(x, 1:20, sep=), 4))
 summary(lm(y ~ x))$r.squared
 [1] 0.2555885
 rsq - function(x, y) {
 + sse - sum(ave(y, x, FUN=function(x) x - mean(x))^2)
 + sstotal - var(y) * (length(y) - 1)
 + 1 - sse / sstotal
 + }
 rsq(x, y)
 [1] 0.2555885
 set.seed(1)
 y - rnorm(8e5)
 x - factor(rep(paste(x, 1:2e4, sep=), 40))
 system.time(rsq(x, y))
 [1] 1.99 0.03 2.06   NA   NA

 Andy


 From: ivo welch

 dear R wizards:

 X is factor with 20,000*20=800,000 observations of 20,000 factors.
 I.e., each factor has 20 observations.  y is 800,000 normally
 distributed data points.  I want to see how much R^2 the X
 factors can provide.  Easy, right?

 lm ( y ~ X)
 and
  aov( y ~ X)
 Error: cannot allocate vector of size 3125000 Kb

 is this computationally infeasible?  (I am not an expert, but
 off-hand, I thought this can work as long as the X's are just
 factors = fitted means.)

 help.search(fixed.effects);

 fixed.effects(nlme) Extract Fixed Effects
 fixed.effects.lmList(nlme)  Extract lmList Fixed Effects
 lme(nlme)   Linear Mixed-Effects Models
 lmeStruct(nlme) Linear Mixed-Effects Structure
 nlme(nlme)  Nonlinear Mixed-Effects Models
 nlmeStruct(nlme)Nonlinear Mixed-Effects
 Structure

 Type 'help(FOO, package = PKG)' to inspect entry 'FOO(PKG) TITLE'.

 ok---I want to read the fixed.effects help.  could the help
 system tell me how to inspect these entries?

 help(fixed.effects)
 wrong
 help.search(fixed.effects)
 two entries, as above.  nothing new.
 ?fixed.effects
 wrong

 eventually, it dawned on me that nlme in parens was not a
 function argument, but the name of the package within which
 fixed.effects lives.  Suggestion:  maybe a different notation
 to name packages would be more intuitive in the help system.
 yes, I know it now, but other novices may not.  even a colon
 instead of a () may be more intuitive in this context.

 library(nlme); ?lme
 and then
 lme(y ~ X)
 Error in getGroups.data.frame(dataMix, groups) :
 Invalid formula for groups


 now I have to beg for some help.  ok, blatant free-riding.
 the lme docs tells me it does the Laird and Ware model, but I
 do not know this model.  the only two examples given at the
 end of the lme help file seem to be similar to what I just
 specified.  so, how do I execute a simple fixed-effects
 model?  (later, I may want to add a variable Z that is a
 continuous random variable.)  could someone please give me
 one quick example?  help is, as always, highly appreciated.

 sincerely,

 /ivo welch

 __
 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

 --
 
 Notice:  This e-mail message, together with any attachments,
 contains information of Merck  Co., Inc. (One Merck Drive,
 Whitehouse Station, New Jersey, USA 08889), and/or its
 affiliates (which may be known outside the United States as
 Merck Frosst, Merck Sharp  Dohme or MSD and in Japan, as
 Banyu) that may be confidential, proprietary copyrighted
 and/or legally privileged. It is intended solely for the use
 of the individual or entity named on this message.  If you
 are not the intended recipient, and have received this
 message in error, please notify us immediately by reply
 e-mail and then delete it from your system.