[R] changing order of lattice plots

2009-03-23 Thread Dan Kortschak
Hi,

This is another question relating to my 2 factor figure.

densityplot(~End-Begin | Type * Chromosome, data=Mon, layout=c(5,12),
xlab="Element Length",type="percent", col="grey60",
strip=strip.custom(style=3, bg="grey90", par.strip.text=list(cex=0.5)))

I would like to flip the plot so those at the bottom are at the top and
so on. I have tried using a `index.cond=list(60:1)' (I have 3 classes
for Chromosome and 20 classes for Type) parameter - this approximates
what I want while I sort out the syntax. But I get `Error in
cond.orders(foo) : Invalid value of index.cond' returned.

The equivalent parameter for a 1x3 plot:

histogram(~Length/1000 | Chromosome, data=readlengths, layout=c(1,3),
xlab="Contig Length (kb)", nint=25, type="count", col="grey60",
index.cond=list(3:1), strip=strip.custom(bg="grey90"))

works fine, so I'm not sure where to go from here.

Any help is greatly appreciated.

thanks
Dan

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


Re: [R] changing order of lattice plots

2009-03-23 Thread Dimitris Rizopoulos

try setting the 'as.table' argument to TRUE, e.g.,

densityplot(~End-Begin | Type * Chromosome, data=Mon, layout=c(5,12), 
xlab="Element Length", type="percent", col="grey60",

strip=strip.custom(style=3, bg="grey90", par.strip.text=list(cex=0.5)),
as.table = TRUE)


I hope it helps.

Best,
Dimitris


Dan Kortschak wrote:

Hi,

This is another question relating to my 2 factor figure.

densityplot(~End-Begin | Type * Chromosome, data=Mon, layout=c(5,12),
xlab="Element Length",type="percent", col="grey60",
strip=strip.custom(style=3, bg="grey90", par.strip.text=list(cex=0.5)))

I would like to flip the plot so those at the bottom are at the top and
so on. I have tried using a `index.cond=list(60:1)' (I have 3 classes
for Chromosome and 20 classes for Type) parameter - this approximates
what I want while I sort out the syntax. But I get `Error in
cond.orders(foo) : Invalid value of index.cond' returned.

The equivalent parameter for a 1x3 plot:

histogram(~Length/1000 | Chromosome, data=readlengths, layout=c(1,3),
xlab="Contig Length (kb)", nint=25, type="count", col="grey60",
index.cond=list(3:1), strip=strip.custom(bg="grey90"))

works fine, so I'm not sure where to go from here.

Any help is greatly appreciated.

thanks
Dan

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



--
Dimitris Rizopoulos
Assistant Professor
Department of Biostatistics
Erasmus University Medical Center

Address: PO Box 2040, 3000 CA Rotterdam, the Netherlands
Tel: +31/(0)10/7043478
Fax: +31/(0)10/7043014

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


Re: [R] If statement generates two outputs

2009-03-23 Thread Wacek Kusnierczyk
jimdare wrote:
> Hi, 
>
> I tried to create the following if / else statement but I keep getting the
> error "Error: unexpected '}' in "size="large",center="none")
> }" (I have highlighted the } in bold where the error is occuring).  I can't
> seem to find a reason for this, does anyone know how I can fix it?
>
> Thanks,
> James  
>
> if (nostocks<=3)
>   
>   {data1<-test[,data1stocks];
>   tex1<-latex(data1, file=paste(i$Species[1], "1.tex", sep=""), 
>rowname = NULL, 
>cgroup = c("Fishstock", stocknames,"Total"), 
>n.cgroup = c(1, rep(2,(nostocks+1)), 
>colheads = c("Year", rep(c("Catch", "TACC"),
> nostocks+1)),
>   size="large",center="none")
>
> }else)
>   

the error message tells you that '}' is unexpected here, and that's
because it really is.  you forgot to close one parenthesis up there:

n.cgroup = c(1, rep(2,(nostocks+1)),

you open three, you close two.  close three and it goes away.

vQ

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


Re: [R] changing order of lattice plots

2009-03-23 Thread Dan Kortschak
Perfect! Thanks.

Dan

On Mon, 2009-03-23 at 08:35 +0100, Dimitris Rizopoulos wrote:
> try setting the 'as.table' argument to TRUE, e.g.,
> 
> densityplot(~End-Begin | Type * Chromosome, data=Mon, layout=c(5,12), 
> xlab="Element Length", type="percent", col="grey60",
> strip=strip.custom(style=3, bg="grey90",
> par.strip.text=list(cex=0.5)),
> as.table = TRUE)
> 
> 
> I hope it helps.
> 
> Best,
> Dimitris

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


Re: [R] Help with 'boot'

2009-03-23 Thread Camille LEPOITTEVIN

Hi Paul,
I struggled with the boot package too, and I found this help link very 
useful:

http://www.mayin.org/ajayshah/KB/R/documents/boot.html
"Getting started with the "boot" package in R for bootstrap inference".
I hope you'll find what you need there !
The second argument for the bootstrap test just tells R on what you want 
to bootstrap on. Just try the examples in this tutorial...


Camille

pgseye a écrit :

Hi,

I'm wanting to test for a difference in medians between 2 groups using
resampling methods. I found the boot package, but don't really understand
how to write the 'statistic' function required as the 2nd argument for the
bootstrap test.

Thanks if you can help,

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


Re: [R] lattice multipanel strip placement - with two factors

2009-03-23 Thread baptiste auguie
I'm not sure I understood your problem (can you provide an  
reproducible example?), but perhaps you can try useOuterStrips() in  
the latticeExtra package (the formatting becomes similar to that of  
the ggplot2 package, perhaps another option to consider)



Hope this helps,

baptiste

On 23 Mar 2009, at 06:54, Dan Kortschak wrote:


Hi,

I'm making a multipanel lattice densityplot figure with 2 factors (3  
and

20 classes in each factor) with the following statement (the
type="percent" is there to prevent plotting the actual points which
detract from the figure - is there another way of doing this?):

densityplot(~End-Begin | Type * Chromosome, data=Mon, layout=c(5,12),
xlab="Element Length",type="percent", col="grey60",
strip=strip.custom(style=3, bg="grey90",  
par.strip.text=list(cex=0.5)))


Plotting 60 panels and associated strips on a page leaves the whole
thing pretty tight and so I'd like to move the 3 class factor strips  
to

the left margin of the whole figure.

Like so (pardon the ASCII art):

+-++++++
+ ++++++
| ||||||
+A++++++
+ ++++++
| ||||||
+-++++++
+ ++++++
| ||||||
+B++++++
+ ++++++
| ||||||
+-++++++
+ ++++++
| ||||||
+C++++++
+ ++++++
| ||||||
+-++++++

Is this possible with lattice and if so, what do I need to do to get  
it.

I've tried strip.left=TRUE and that just make the problem a horizontal
one rather than a vertical one.

Any help is greatly appreciated.

thanks
Dan

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


_

Baptiste Auguié

School of Physics
University of Exeter
Stocker Road,
Exeter, Devon,
EX4 4QL, UK

Phone: +44 1392 264187

http://newton.ex.ac.uk/research/emag

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


[R] performance: zoo's rollapply() vs inline

2009-03-23 Thread Ken-JP

zoo's rollapply() function appears to be extremely useful for plugging in a
function on-the-fly to run over a window.  With inline, there is a lot more
coding and room for error, and the code is less portable because the user
has to have R compiling set up or it won't work.

However, rollapply() seems to be really slow.  Several orders of magnitude
slower than inline, in fact.  I don't know how to call R functions from C
inline yet, but it looks like I need to learn, because the speed difference
is just way too big.

The results of a quick test are shown below.  

I am totally open to suggestions on how to do windowed calculations, in
general, but it looks like I may have to bite the bullet and learn all the
intricacies of calling R from C.

NOTE:  pchg.inline() is not shown because it's much longer/complex than
pchg.rollapply(), but I am doing no optimizations.



pchg.rollapply <-   function(this, m, shift=1, ...)  {
  rollapply( m, shift+1, function(x) { x[shift+1]/x[1] - 1; }, align="right"
);
}

> dim( m )
[1] 4518  800
> system.time( x.rollapply <- pchg.rollapply( m, 20 ) )
   user  system elapsed 
 146.940.81  157.03 
> system.time( x.inline<- pchg.inline( m, 20 ) )
   user  system elapsed 
   0.690.000.72 



-- 
View this message in context: 
http://www.nabble.com/performance%3A--zoo%27s-rollapply%28%29-vs-inline-tp22656214p22656214.html
Sent from the R help mailing list archive at Nabble.com.

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


[R] Syntax in lmer...

2009-03-23 Thread Anne Kempel
Sorry for bothering you - but I have a little problem. I`m completely 
new in R and trying to use lmer. That works, but I m not sure if my code 
is right, so I would just like to get an opinion from the experts.


I have a Growthrate of plants (RGR
Fixed Factors: Fertilizer, Status (alien or native plants) and their 
interaction
Random Factors: Age is a covariable; than I have plantfamily and Species 
nested within Family.
Now I want also to test for the interaction of Fertilizer with 
plantfamily and species. Did I do this the right way, and is the 0+ 
necessary?


model<-lmer(RGR~Fertilizer+Status+Fertilizer:Status+(1|Age)+(1|Family/Species)+(0+Fertilizer|Family/Species),data=rgr)

I`m thankful for any comments!
anne



--
Anne Kempel
Institute of Plant Sciences
University of Bern
Altenbergrain 21
CH-3103 Bern
kem...@ips.unibe.ch

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


[R] How to get commands history as a character vector instead of displaying them?

2009-03-23 Thread Yihui Xie
Hi Everyone,

I want to get the commands history as a character vector instead of
just displaying them, but the function history() just returns NULL. I
checked the source code of 'history' and could not find a solution.
Anybody has an idea? Thanks!

P. S. My original problem is, when a user opens a graphics device like
png() or pdf(), I want to know the file name used by this device. I
thought history() would help, but it could not.

> sessionInfo()
R version 2.8.1 (2008-12-22)
i386-pc-mingw32

locale:
LC_COLLATE=Chinese_People's Republic of
China.936;LC_CTYPE=Chinese_People's Republic of
China.936;LC_MONETARY=Chinese_People's Republic of
China.936;LC_NUMERIC=C;LC_TIME=Chinese_People's Republic of China.936

attached base packages:
[1] stats graphics  grDevices utils datasets  methods   base

Regards,
Yihui
--
Yihui Xie 
Phone: +86-(0)10-82509086 Fax: +86-(0)10-82509086
Mobile: +86-15810805877
Homepage: http://www.yihui.name
School of Statistics, Room 1037, Mingde Main Building,
Renmin University of China, Beijing, 100872, China

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


Re: [R] performance: zoo's rollapply() vs inline

2009-03-23 Thread Achim Zeileis

On Mon, 23 Mar 2009, Ken-JP wrote:


zoo's rollapply() function appears to be extremely useful for plugging in a
function on-the-fly to run over a window.  With inline, there is a lot more
coding and room for error, and the code is less portable because the user
has to have R compiling set up or it won't work.

However, rollapply() seems to be really slow.  Several orders of magnitude
slower than inline, in fact.  I don't know how to call R functions from C
inline yet, but it looks like I need to learn, because the speed difference
is just way too big.


It depends what you want to do with it. If you use rollapply() for 
operations that you could do in a vectorized way then it is certainly not 
a good idea (see below). Important functions, especially rolling means, 
are special cased and are much faster than a regular rollapply().



The results of a quick test are shown below.

I am totally open to suggestions on how to do windowed calculations, in
general, but it looks like I may have to bite the bullet and learn all the
intricacies of calling R from C.

NOTE:  pchg.inline() is not shown because it's much longer/complex than
pchg.rollapply(), but I am doing no optimizations.



pchg.rollapply <-   function(this, m, shift=1, ...)  {
 rollapply( m, shift+1, function(x) { x[shift+1]/x[1] - 1; }, align="right"
);
}


This is really a bad example because your function is flawed (no 
dependence on "this") and it is not clear to me why you would want to use 
rollapply(). Just doing

  m/lag(m, -shift) - 1
should do the job.
Z


dim( m )

[1] 4518  800

system.time( x.rollapply <- pchg.rollapply( m, 20 ) )

  user  system elapsed
146.940.81  157.03

system.time( x.inline<- pchg.inline( m, 20 ) )

  user  system elapsed
  0.690.000.72



--
View this message in context: 
http://www.nabble.com/performance%3A--zoo%27s-rollapply%28%29-vs-inline-tp22656214p22656214.html
Sent from the R help mailing list archive at Nabble.com.

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




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


Re: [R] How to get commands history as a character vector instead of displaying them?

2009-03-23 Thread Romain Francois

Yihui Xie wrote:

Hi Everyone,

I want to get the commands history as a character vector instead of
just displaying them, but the function history() just returns NULL. I
checked the source code of 'history' and could not find a solution.
Anybody has an idea? Thanks!
  
history eventually calls file.show, which will use the pager option to 
determine how to show the file, so you can do something like that:


history <- function( ... ){
old.op <- options( pager = function( files, header, title, delete.file ) 
readLines( files ) ); on.exit( options( old.op ) )

utils::history(...)
}
history( pattern = "png" )

I don't see a way to get device information other than the name of the 
device (with dev.cur)


Romain


P. S. My original problem is, when a user opens a graphics device like
png() or pdf(), I want to know the file name used by this device. I
thought history() would help, but it could not.

  

sessionInfo()


R version 2.8.1 (2008-12-22)
i386-pc-mingw32

locale:
LC_COLLATE=Chinese_People's Republic of
China.936;LC_CTYPE=Chinese_People's Republic of
China.936;LC_MONETARY=Chinese_People's Republic of
China.936;LC_NUMERIC=C;LC_TIME=Chinese_People's Republic of China.936

attached base packages:
[1] stats graphics  grDevices utils datasets  methods   base

Regards,
Yihui
--
Yihui Xie 
Phone: +86-(0)10-82509086 Fax: +86-(0)10-82509086
Mobile: +86-15810805877
Homepage: http://www.yihui.name
School of Statistics, Room 1037, Mingde Main Building,
Renmin University of China, Beijing, 100872, China

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


  



--
Romain Francois
Independent R Consultant
+33(0) 6 28 91 30 30
http://romainfrancois.blog.free.fr

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


Re: [R] Syntax in lmer...

2009-03-23 Thread ONKELINX, Thierry
Dear Anne,

Is age a continuous or a categorical variable? If it is a continuous
variable, then I would add it to the fixed effects.

You have two options for the random effect: (Fertilier|Family/Species)
or (1|Family/Species) + (0 + Fertilizer|Family/Species). Both yield a
random intercept and a random slope along Fertilizer. But the main
difference is the correlation between the random effects. The first
model allows for correlation between the random effect and the random
slope. Whereas the second model implies that they are independent.

This is what I would do to test for the interaction.
#If age is categorical
M0 <- lmer(RGR ~ Fertilizer * Status + (1|Age) + (1|Family/Species) + (0
+ Fertilizer|Family/Species), data=rgr, REML =  TRUE)
M1 <- lmer(RGR ~ Fertilizer * Status + (1|Age) + (1|Family/Species),
data=rgr, REML =  TRUE)
anova(M0, M1)

#If age is continuous
M0 <- lmer(RGR ~ Fertilizer * Status + Age + (1|Family/Species) + (0 +
Fertilizer|Family/Species), data=rgr, REML =  TRUE)
M1 <- lmer(RGR ~ Fertilizer * Status + Age + (1|Family/Species),
data=rgr, REML =  TRUE)
anova(M0, M1)

HTH,

Thierry

PS. The R-Sig-mixed models mailing list is better suited for questions
on lme4.



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

To call in the statistician after the experiment is done may be no more
than asking him to perform a post-mortem examination: he may be able to
say what the experiment died of.
~ Sir Ronald Aylmer Fisher

The plural of anecdote is not data.
~ Roger Brinner

The combination of some data and an aching desire for an answer does not
ensure that a reasonable answer can be extracted from a given body of
data.
~ John Tukey

-Oorspronkelijk bericht-
Van: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org]
Namens Anne Kempel
Verzonden: maandag 23 maart 2009 9:54
Aan: r-help@r-project.org
Onderwerp: [R] Syntax in lmer...

Sorry for bothering you - but I have a little problem. I`m completely 
new in R and trying to use lmer. That works, but I m not sure if my code

is right, so I would just like to get an opinion from the experts.

I have a Growthrate of plants (RGR
Fixed Factors: Fertilizer, Status (alien or native plants) and their 
interaction
Random Factors: Age is a covariable; than I have plantfamily and Species

nested within Family.
Now I want also to test for the interaction of Fertilizer with 
plantfamily and species. Did I do this the right way, and is the 0+ 
necessary?

model<-lmer(RGR~Fertilizer+Status+Fertilizer:Status+(1|Age)+(1|Family/Sp
ecies)+(0+Fertilizer|Family/Species),data=rgr)

I`m thankful for any comments!
anne



-- 
Anne Kempel
Institute of Plant Sciences
University of Bern
Altenbergrain 21
CH-3103 Bern
kem...@ips.unibe.ch

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

Dit bericht en eventuele bijlagen geven enkel de visie van de schrijver weer 
en binden het INBO onder geen enkel beding, zolang dit bericht niet bevestigd is
door een geldig ondertekend document. The views expressed in  this message 
and any annex are purely those of the writer and may not be regarded as stating 
an official position of INBO, as long as the message is not confirmed by a duly 
signed document.

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


Re: [R] lattice multipanel strip placement - with two factors

2009-03-23 Thread Dan Kortschak
That's part way to my desired solution, but not completely there.

Here's an example:

chromosome<-rep(c("A","X","Y"),time=20)
type<-rep(c(1:20),times=3)
length<-rnorm(60)
densityplot(~length | type * chromosome, layout=c(5,12))

What I would like to see is the chromosome strip (A, X, Y) once on the
left as useOuterStrips() gives, but with the type strip internal as the
normal lattice function gives.

using 

useOuterStrips(densityplot(~length | type * chromosome, layout=c(5,12)))

for my data (data ranges from 0-400) is about as unreadable at the
normal situation (the toy example doesn't get that across).

thanks
Dan

On Mon, 2009-03-23 at 08:40 +, baptiste auguie wrote:
> I'm not sure I understood your problem (can you provide an  
> reproducible example?), but perhaps you can try useOuterStrips() in  
> the latticeExtra package (the formatting becomes similar to that of  
> the ggplot2 package, perhaps another option to consider)

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


Re: [R] How to get commands history as a character vector instead of displaying them?

2009-03-23 Thread Erich Neuwirth
Your original question probably is answered by
dev.cur and dev.list




Yihui Xie wrote:
> Hi Everyone,
> 
> I want to get the commands history as a character vector instead of
> just displaying them, but the function history() just returns NULL. I
> checked the source code of 'history' and could not find a solution.
> Anybody has an idea? Thanks!
> 
> P. S. My original problem is, when a user opens a graphics device like
> png() or pdf(), I want to know the file name used by this device. I
> thought history() would help, but it could not.
> 
>> sessionInfo()
> R version 2.8.1 (2008-12-22)
> i386-pc-mingw32
> 
> locale:
> LC_COLLATE=Chinese_People's Republic of
> China.936;LC_CTYPE=Chinese_People's Republic of
> China.936;LC_MONETARY=Chinese_People's Republic of
> China.936;LC_NUMERIC=C;LC_TIME=Chinese_People's Republic of China.936
> 
> attached base packages:
> [1] stats graphics  grDevices utils datasets  methods   base
> 
> Regards,
> Yihui
> --
> Yihui Xie 
> Phone: +86-(0)10-82509086 Fax: +86-(0)10-82509086
> Mobile: +86-15810805877
> Homepage: http://www.yihui.name
> School of Statistics, Room 1037, Mingde Main Building,
> Renmin University of China, Beijing, 100872, China
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 
> 
> 
> 
> 
> No virus found in this incoming message.
> Checked by AVG - www.avg.com 
> Version: 8.0.238 / Virus Database: 270.11.24/2017 - Release Date: 03/22/09 
> 17:51:00
> 

-- 
Erich Neuwirth, University of Vienna
Faculty of Computer Science
Computer Supported Didactics Working Group
Visit our SunSITE at http://sunsite.univie.ac.at
Phone: +43-1-4277-39464 Fax: +43-1-4277-39459

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


Re: [R] performance: zoo's rollapply() vs inline

2009-03-23 Thread Gabor Grothendieck
On Mon, Mar 23, 2009 at 5:08 AM, Achim Zeileis
 wrote:
> On Mon, 23 Mar 2009, Ken-JP wrote:
>
>> zoo's rollapply() function appears to be extremely useful for plugging in
>> a
>> function on-the-fly to run over a window.  With inline, there is a lot
>> more
>> coding and room for error, and the code is less portable because the user
>> has to have R compiling set up or it won't work.
>>
>> However, rollapply() seems to be really slow.  Several orders of magnitude
>> slower than inline, in fact.  I don't know how to call R functions from C
>> inline yet, but it looks like I need to learn, because the speed
>> difference
>> is just way too big.
>
> It depends what you want to do with it. If you use rollapply() for
> operations that you could do in a vectorized way then it is certainly not a
> good idea (see below). Important functions, especially rolling means, are
> special cased and are much faster than a regular rollapply().
>
>> The results of a quick test are shown below.
>>
>> I am totally open to suggestions on how to do windowed calculations, in
>> general, but it looks like I may have to bite the bullet and learn all the
>> intricacies of calling R from C.
>>
>> NOTE:  pchg.inline() is not shown because it's much longer/complex than
>> pchg.rollapply(), but I am doing no optimizations.
>>
>>
>> 
>>
>> pchg.rollapply <-   function(this, m, shift=1, ...)  {
>>  rollapply( m, shift+1, function(x) { x[shift+1]/x[1] - 1; },
>> align="right"
>> );
>> }
>
> This is really a bad example because your function is flawed (no dependence
> on "this") and it is not clear to me why you would want to use rollapply().
> Just doing
>  m/lag(m, -shift) - 1
> should do the job.

Also check out ?diff.zoo and note the arithmetic = FALSE argument although
using lag.zoo as above may be slightly faster.

If the above is not your real problem also check out rollmean, rollmedian
and rollmax which are also written in R but have been optimized for those
operations.

The xts package has some routines written in C that may also be applicable.

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


[R] Efficiently create dummy

2009-03-23 Thread Rob Denniker
What's the neat way to create a dummy from a list?
The code below is not replicable, but hopefully self-explanatory...

d$treatment<-rep(1,length(d))

notreat<-c("AR", "DE", "MS", "NY", "TN", "AK", "LA", "MD",  "NC", "OK", "UT", 
"VA")

#i would really like this to work:
d$treatment[d$st==any(notreat)]<-0

#but instead i resort to this
for (i in 1:length(notreat)) {
temp.st <- notreat[i]
d$treatment[d$st==temp.st]<-0
i<-i+1 }


Thanks, list!

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


Re: [R] Efficiently create dummy

2009-03-23 Thread Dimitris Rizopoulos

well, you can first define a factor, and the use model.matrix(), e.g.,

fc <- factor(c("AR", "DE", "MS", "NY"))
model.matrix(~ fc)


for more info check ?model.matrix(), e.g., if you want to change the 
contrasts.


I hope it helps.

Best,
Dimitris


Rob Denniker wrote:

What's the neat way to create a dummy from a list?
The code below is not replicable, but hopefully self-explanatory...

d$treatment<-rep(1,length(d))

notreat<-c("AR", "DE", "MS", "NY", "TN", "AK", "LA", "MD",  "NC", "OK", "UT", 
"VA")

#i would really like this to work:
d$treatment[d$st==any(notreat)]<-0

#but instead i resort to this
for (i in 1:length(notreat)) {
temp.st <- notreat[i]
d$treatment[d$st==temp.st]<-0
i<-i+1 }


Thanks, list!

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



--
Dimitris Rizopoulos
Assistant Professor
Department of Biostatistics
Erasmus University Medical Center

Address: PO Box 2040, 3000 CA Rotterdam, the Netherlands
Tel: +31/(0)10/7043478
Fax: +31/(0)10/7043014

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


Re: [R] xlsReadWrite library

2009-03-23 Thread Hans-Peter Suter
2009/3/19 Pascal Candolfi :
> Any idea why this library was removed and where could I find it for Windows
> (only Unix in the Archive) ?

As indicated it was for for the binary (non FOSS) component.

Yesterday I checked the existing packages (xlsReadWrite(Pro) on
http://treetron.googlepages.com) if they run on 'R 2.9.0 under
development'. They do fine (well apart from two or three known _small_
bugs ;-).

I made a fully OSS and C-based version meanwhile which downloads (or
gives instruction how to compile for yourself) the necessary DLL from
a non-CRAN source.

However as I wanted to make sure that this 'packaging' still counts as
a 'binary distribution' I explicitely asked the developer from the 3rd
library (Flexcel) if this is ok. Last Friday I got his approval, so at
least technically/legally it is possibled.

When I'm happy with the new RUnit tests, adapted build scripts, etc.,
I'll ask again if 'xlsReadWrite' will be accepted on CRAN in this new
source-only form. I hope so, but if not then there are other
possibilities. In any case I intend to put the repo at some public
space (most probably github).

Cheers,
Hans-Peter

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


Re: [R] Efficiently create dummy

2009-03-23 Thread Chuck Cleland
On 3/23/2009 5:44 AM, Rob Denniker wrote:
> What's the neat way to create a dummy from a list?
> The code below is not replicable, but hopefully self-explanatory...
> 
> d$treatment<-rep(1,length(d))
> 
> notreat<-c("AR", "DE", "MS", "NY", "TN", "AK", "LA", "MD",  "NC", "OK", "UT", 
> "VA")
> 
> #i would really like this to work:
> d$treatment[d$st==any(notreat)]<-0
> 
> #but instead i resort to this
> for (i in 1:length(notreat)) {
> temp.st <- notreat[i]
> d$treatment[d$st==temp.st]<-0
> i<-i+1 }

d$treatment <- as.numeric(!(d$st %in% notreat))

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

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

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


Re: [R] How to get commands history as a character vector instead of displaying them?

2009-03-23 Thread Wacek Kusnierczyk
Romain Francois wrote:
> Yihui Xie wrote:
>> Hi Everyone,
>>
>> I want to get the commands history as a character vector instead of
>> just displaying them, but the function history() just returns NULL. I
>> checked the source code of 'history' and could not find a solution.
>> Anybody has an idea? Thanks!
>>   
> history eventually calls file.show, which will use the pager option to
> determine how to show the file, so you can do something like that:
>
> history <- function( ... ){
> old.op <- options( pager = function( files, header, title, delete.file
> ) readLines( files ) ); on.exit( options( old.op ) )
> utils::history(...)
> }
> history( pattern = "png" )

i think the following is an acceptable alternative:

history = function() {
   file = tempfile()
   on.exit(unlink(file))
   savehistory(file)
   readLines(file) }

the output is *lines* of text, but if you need whole executable
expressions, you can parse the output:

1 + 1
ph = parse(text=history())
as.list(ph)
ph[length(ph)-1]
# expression(1 + 1)
eval(ph[length(ph)-1])
# [1] 2


vQ

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


Re: [R] fisher.test - FEXACT error 7

2009-03-23 Thread Bernardo Rangel Tura
On Fri, 2009-03-20 at 18:29 +, Helena Mouriño wrote:
> Dear all,
> 
> Im having an awkward problem in R.  When I write the command
> Fisher.test(school.data,workspace=2e+07), where school.data is the matrix
> corresponding to the data set,
> I get the error  message:

> FEXACT error 7.
> LDSTP is too small for this problem.
> Try increasing the size of the workspace.
> 
> Increasing the workspace: 
> Fisher.test(school.data,workspace=1e+10),
> I get a different message, but it still doesnt work:
> 
> NAs in foreign function call (arg 10)
> In addition: Warning message:
> In fisher.test(dados, workspace = 1e+10, alternative = "two.sided") :
>   NAs introduced by coercion

Hi Helena,

In this case you can try 3 solutions:

1- chisq.test(school.data), but pay attention if expected value of any
cell is < 5 

2- Fisher.test(school.data,workspace=2e+07,hybrid=TRUE) from Help

For larger than  2 by 2 tables and 'hybrid = TRUE', asymptotic
 chi-squared probabilities are only used if the 'Cochran
 conditions' are satisfied, that is if no cell has count zero, and
 more than 80% of the cells have counts at least 5.

3- Use "large tables" approach from Sir David Cox:

Law, G. R. and Cox, D. R. and Machonochie, N. E. S. and Simpson, J. and
Roman, E. and Carpenter, L. M. (2001) Large tables. Biostatistics
2(2):pp. 163-171.


-- 
Bernardo Rangel Tura, M.D,MPH,Ph.D
National Institute of Cardiology
Brazil

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


[R] subset and "or" operator

2009-03-23 Thread lauramorg...@bluewin.ch
Hello,
I'm trying to subset a dataframe where I have many observation taken in 
different years.
I would like to subset the dataframe (in this example called "table") to get a 
new dataframe containing only the 
observation of year 1995, 1998 and 2000.
I've tried to use subset and the or operator "|" like this:
subset(table, year==1995|1998|2000)-->table2
but the dataframe didn't subset and I get an object called table2 that is 
identical to table.
I guess I'm doing something wrong!!

Thank in advance for any suggestions!
Laura

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


Re: [R] any package for connecting berkeley db in R?

2009-03-23 Thread Bernardo Rangel Tura
On Fri, 2009-03-20 at 17:23 -0400, Zheng, Xin (NIH) [C] wrote:
> Hi there,
> 
> Is there any such package? I searched but found none. Thanks in advance.
> 
> Xin Zheng
> 

Hi Xin

Do you try package DBI?

-- 
Bernardo Rangel Tura, M.D,MPH,Ph.D
National Institute of Cardiology
Brazil

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


Re: [R] newton method

2009-03-23 Thread Rau, Roland
Hi,

you might be also interested in a general overview as given here:
http://cran.r-project.org/web/views/Optimization.html

Hope this helps,
Roland


> -Original Message-
> From: r-help-boun...@r-project.org 
> [mailto:r-help-boun...@r-project.org] On Behalf Of Michael Kubovy
> Sent: Monday, March 23, 2009 4:35 AM
> To: Roslina Zakaria
> Cc: r-help@r-project.org
> Subject: Re: [R] newton method
> 
> Take a look at the functionsnlm(), optim() in the stats package and  
> maxNR() in the maxLik package.
> 
> On Mar 22, 2009, at 11:15 PM, Roslina Zakaria wrote:
> 
> > Does R has a topic on newton's method?
> 
> 
> _
> Professor Michael Kubovy
> University of Virginia
> Department of Psychology
> Postal Address:
>   P.O.Box 400400, Charlottesville, VA 22904-4400
> Express Parcels Address:
>   Gilmer Hall, Room 102, McCormick Road, Charlottesville, VA 22903
> Office:B011;  Phone: +1-434-982-4729
> Lab:B019; Phone: +1-434-982-4751
> WWW:http://www.people.virginia.edu/~mk9y/
> Skype name: polyurinsane
> 
> 
> 
> 
> 
>   [[alternative HTML version deleted]]
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide 
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 

--
This mail has been sent through the MPI for Demographic Research.  Should you 
receive a mail that is apparently from a MPI user without this text displayed, 
then the address has most likely been faked. If you are uncertain about the 
validity of this message, please check the mail header or ask your system 
administrator for assistance.

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


Re: [R] subset and "or" operator

2009-03-23 Thread Erich Neuwirth
subset(table, year %in% c(1995,1998,2000))-->table2

> Hello,
> I'm trying to subset a dataframe where I have many observation taken in 
> different years.
> I would like to subset the dataframe (in this example called "table") to get 
> a new dataframe containing only the 
> observation of year 1995, 1998 and 2000.
> I've tried to use subset and the or operator "|" like this:
> subset(table, year==1995|1998|2000)-->table2
> but the dataframe didn't subset and I get an object called table2 that is 
> identical to table.
> I guess I'm doing something wrong!!
> 
> Thank in advance for any suggestions!
> Laura
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 
> 
> 
> 
> 
> No virus found in this incoming message.
> Checked by AVG - www.avg.com 
> Version: 8.0.238 / Virus Database: 270.11.24/2018 - Release Date: 03/23/09 
> 06:52:00
> 

-- 
Erich Neuwirth, University of Vienna
Faculty of Computer Science
Computer Supported Didactics Working Group
Visit our SunSITE at http://sunsite.univie.ac.at
Phone: +43-1-4277-39464 Fax: +43-1-4277-39459

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


[R] Capitalizing first letter of word or phrase

2009-03-23 Thread Daren Tan
I managed to find toupper() which translates all letters to uppercase.
Is there a function to capitalize only the first letter of word or
phrase ?

Thanks

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


Re: [R] Capitalizing first letter of word or phrase

2009-03-23 Thread Gabor Grothendieck
See examples on ?toupper page.

On Mon, Mar 23, 2009 at 8:03 AM, Daren Tan  wrote:
> I managed to find toupper() which translates all letters to uppercase.
> Is there a function to capitalize only the first letter of word or
> phrase ?
>
> Thanks
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

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


[R] Iterative Proportional Fitting, use

2009-03-23 Thread Koen Hufkens
Hi list,

I would like to normalize a matrix (two actually for comparison) using
iterative proportional fitting.

Using ipf() would be the easiest way to do this, however I can't get my
head around the use of the function. More specifically, the margins
settings...

for a matrix:

mat <- matrix(c(65,4,22,24,6,81,5,8,0,11,85,19,4,7,3,90),4,4)

using

fit <- ipf(mat,margins=c(1,1,1,1,0,1,1,1,1))

generates a matrix with just 1's.

using


fit <- ipf(mat,margins=c(100,100,100,100,0,100,100,100,100))

gives a segmentation fault and crashes R !

so how do you define the margin values to which to sum the row and
column values in your matrix correctly?

Kind regards,
Koen

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


Re: [R] fisher.test - FEXACT error 7

2009-03-23 Thread Frank E Harrell Jr

Bernardo Rangel Tura wrote:

On Fri, 2009-03-20 at 18:29 +, Helena Mouriño wrote:

Dear all,

Im having an awkward problem in R.  When I write the command
Fisher.test(school.data,workspace=2e+07), where school.data is the matrix
corresponding to the data set,
I get the error  message:



FEXACT error 7.
LDSTP is too small for this problem.
Try increasing the size of the workspace.

Increasing the workspace: 
Fisher.test(school.data,workspace=1e+10),

I get a different message, but it still doesnt work:

NAs in foreign function call (arg 10)
In addition: Warning message:
In fisher.test(dados, workspace = 1e+10, alternative = "two.sided") :
  NAs introduced by coercion


Hi Helena,

In this case you can try 3 solutions:

1- chisq.test(school.data), but pay attention if expected value of any
cell is < 5 


Note that this requirement was never really checked by Pearson (as 
pointed out by Cochran) and is overly cautious.  Generally the 
Pearson-Cochran chi-square test works well much more frequently than 
previously thought, and usually better than Fisher's exact test.  See 
the reference below. -Frank


@Article{cam07chi,
  author =   {Campbell, Ian},
  title ={Chi-squared and {Fisher-Irwin} tests of 
two-by-two tab

les with small sample recommendations},
  journal =  Stat in Med,
  year = 2007,
  volume =   26,
  pages ={3661-3675},
  annote =   {2x2 table;chi-squared test;Fisher-Irwin 
test;exact tests;small sample recommendations;latest edition of 
Armitage's book recommends that continuity adjustments never be used for 
contingency table chi-square tests;E. Pearson modification of Pearson 
chi-square test, differing from the original by a factor of 
(N-1)/N;Cochran noted that the number 5 in "expected frequency less than 
5" was arbitrary;findings of published studies may be summarized as 
follows, for comparative trials:``1. Yate's chi-squared test has type I 
error rates less than the nominal, often less than half the nominal; 2. 
The Fisher-Irwin test has type I error rates less than the nominal; 3. K 
Pearson's version of the chi-squared test has type I error rates closer 
to the nominal than Yate's chi-squared test and the Fisher-Irwin test, 
but in some situations gives type I errors appreciably larger than the 
nominal value; 4. The 'N-1' chi-squared test, behaves like K. Pearson's 
'N' version, but the tendency for higher than nominal values is reduced; 
5. The two-sided Fisher-Irwin test using Irwin's rule is less 
conservative than the method doubling the one-sided probability; 6. The 
mid-P Fisher-Irwin test by doubling the one-sided probability performs 
better than standard versions of the Fisher-Irwin test, and the mid-P 
method by Irwin's rule performs better still in having actual type I 
errors closer to nominal levels."; strong support for the 'N-1' test 
provided expected frequencies exceed 1;flaw in Fisher test which was 
based on Fisher's premise that marginal totals carry no useful 
information;demonstration of their useful information in very small 
sample sizes;Yate's continuity adjustment of N/2 is a large over 
correction and is inappropriate;counter arguments exist to the use of 
randomization tests in randomized trials;calculations of worst 
cases;overall recommendation: use the 'N-1' chi-square test when all 
expected frequencies are at least 1, otherwise use the Fisher-Irwin test 
using Irwin's rule for two-sided tests, taking tables from either tail 
as likely, or less, as that observed; see letter to the editor by 
Antonio Andres and author's reply in 27:1791-1796; 2008.}

}



2- Fisher.test(school.data,workspace=2e+07,hybrid=TRUE) from Help

For larger than  2 by 2 tables and 'hybrid = TRUE', asymptotic
 chi-squared probabilities are only used if the 'Cochran
 conditions' are satisfied, that is if no cell has count zero, and
 more than 80% of the cells have counts at least 5.

3- Use "large tables" approach from Sir David Cox:

Law, G. R. and Cox, D. R. and Machonochie, N. E. S. and Simpson, J. and
Roman, E. and Carpenter, L. M. (2001) Large tables. Biostatistics
2(2):pp. 163-171.





--
Frank E Harrell Jr   Professor and Chair   School of Medicine
 Department of Biostatistics   Vanderbilt University

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


Re: [R] Problems with combining plots

2009-03-23 Thread johnhj

I have still the same problem... As you said I tried with par(mfrow=c(2,1))
and par(mfrow=c(1,2)) but without success. Could R compiler be the problem ?



Wills, Kellie wrote:
> 
> par(mfrow=c(1,1)) will give you just one panel.  Try par(mfrow=c(2,1)) or
> par(mfrow=c(1,2)).
> 
> 
> -Original Message-
> From: r-help-boun...@r-project.org on behalf of johnhj
> Sent: Sun 3/22/2009 10:50 AM
> To: r-help@r-project.org
> Subject: [R]  Problems with combining plots
>  
> 
> Hii,
> 
> I will combine some plots. Like this example here
> http://www.statmethods.net/advgraphs/layout.html I tired to do it for 2
> plots but without success.
> 
> Here is my code:
> 
> test<-read.table(file="D:/file.txt")
> space<-read.table(file="D:/space.txt")
> 
> space$gruppe <- 502*rep(1:6, each=7)
> x<- c(test$V1)
> y<- c(test$V2)
> 
> par(mfrow=c(1,1)) 
> 
> png(filename = "D:/example.png", width = 640, height = 480, pointsize =
> 12,
> bg = "white",  res = NA)
> 
> boxplot(V2 ~ gruppe , data = space , col = "lightgray",boxwex=0.2)
> plot(panel.first=grid(ny=NULL,nx=NULL),x,y, xlab = "Zeit(sec)", ylab
> ="Datenrate(MBit(sec))",ylim=c(0,40), col ="purple", type ="l", main
> ="combined plots",lwd=2)
> 
> dev.off()
> 
> 
> What is the mistake in my code ? 
> 
> 
> -- 
> View this message in context:
> http://www.nabble.com/Problems-with-combining-plots-tp22646692p22646692.html
> Sent from the R help mailing list archive at Nabble.com.
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 
> 

-- 
View this message in context: 
http://www.nabble.com/Problems-with-combining-plots-tp22646692p22658673.html
Sent from the R help mailing list archive at Nabble.com.

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


[R] help

2009-03-23 Thread Dominique Katshunga
Dear R-users,
Can someone help understand why in the following example R would return a 
negative value for the last parameter when I have constrained it to be positive 
and greater than or equal to 0.01?
 
> nlminb(start=initial value,objective=some function, 
> lower=c(-Inf,-Inf,-Inf,0.01),upper=c(Inf,Inf,Inf,Inf))
Many thanks,
Dominique
 
 
Dominique Katshunga
==
Lecturer in Statistical Sciences Department
University of Cape Town
Tel. 021 6504669
Fax. 021 6504773
email. dominique.katshu...@uct.ac.za

[[alternative HTML version deleted]]

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


Re: [R] Iterative Proportional Fitting, use

2009-03-23 Thread Gerard M. Keogh
Keon,

why not fit a loglinear independence model which as far as I know is the
same.

Gerard

Here's an example from Agresti - Intro to Cat Data analysis
Example: Alcohol, cigarette, marijuana use
|--+--+|
|  Alcohol | Cigarette|Marijuana Use   |
|  |  ||
|use   |use   ||
|--+--+|
|  |  |   Yes No   |
|--+--+|
|Yes   |Yes   |   911 538  |
|--+--+|
|  |No|   44 456   |
|--+--+|
|No|Yes   |3 43|
|--+--+|
|  |No|2 279   |
|--+--+|



Coding and Models

  table8.3 = read.table(textConnection("alc cig mar count
  Yes Yes Yes 911
  Yes Yes No  538
  Yes No  Yes 44
  Yes No  No  456
  No  Yes Yes 3
  No  Yes No  43
  No  No  Yes 2
  No  No  No  279"),header=TRUE)
  closeAllConnections()
  # independence model (A,C,M)
  fit1.a.c.m = glm(count ~ mar+cig+alc, family=poisson, data=table8.3)
  fit1.glm$fitted.values
  # intermediate model
  fit2.m.ca = glm(count ~ mar+cig:alc, family=poisson, data=table8.3)
  fit2.m.ca$fitted.values
  # homogeneous association model
  fit3.m.c.a  =  glm(count  ~  mar:cig+mar:alc+cig:alc, family=poisson,
  data=table8.3)
  fit3.m.c.a$fitted.values
  # saturated model
  fits = glm(count ~ mar*cig*alc, family=poisson, data=table8.3)
  fits$fitted.values


The  coding  for  variables  in the above program and the fitted values are

given  below – they show that the homogeneous association model is the only

model that fits these data well.

|-+++++-+---|
| Alcohol |  Cigarette |  Marijuana | Actual |   (A,C,M)  |  (AC,M) | 
(AC:AM:CM)|
|   use   | Use| Use|  (ACM) | Independenc| | 
homogeneou|
| ||||  e | | s 
|
|-+++++-+---|
|   Yes   | Yes| Yes|   911  |540.0   |  611.2  |   910.4   
|
|-+++++-+---|
| || No |   538  |740.2   |  837.8  |   538.6   
|
|-+++++-+---|
| | No | Yes|   44   |282.1   |  210.9  |44.6   
|
|-+++++-+---|
| || No |   456  |386.7   |  289.1  |   455.4   
|
|-+++++-+---|
|No   | Yes| Yes|3   |90.6|   19.4  |3.6
|
|-+++++-+---|
| || No |   43   |124.2   |   26.6  |42.4   
|
|-+++++-+---|
| | No | Yes|2   |47.3|  118.5  |1.4
|
|-+++++-+---|
| || No |   279  |64.9|  162.5  |   279.6   
|
|-+++++-+---|







   
 Koen Hufkens  
  To 
 Sent by:  r-help
 r-help-boun...@r-  cc 
 project.org   
   Subject 
   [R] Iterative Proportional Fitting, 
 23/03/2009 12:13  use 
   
   
 

Re: [R] Iterative Proportional Fitting, use

2009-03-23 Thread Koen Hufkens
The data I used was just an example to work upon.

My real dataset is a confusion matrix of 24x24 (and 17x17), so coding it
into a model with all different kinds of combinations seems tedious.
That's why I hoped to use the ipf() function as it accepts a matrix as
input.

Thanks for the suggestion but I hope I can get around this without
recoding the original data.

Kind regards,
Koen

-Original Message-
From: Gerard M. Keogh 
To: Koen Hufkens 
Cc: r-help , r-help-boun...@r-project.org
Subject: Re: [R] Iterative Proportional Fitting, use
Date: Mon, 23 Mar 2009 13:11:15 +

Keon,

why not fit a loglinear independence model which as far as I know is the
same.

Gerard

Here's an example from Agresti - Intro to Cat Data analysis
Example: Alcohol, cigarette, marijuana use
|--+--+|
|  Alcohol | Cigarette|Marijuana Use   |
|  |  ||
|use   |use   ||
|--+--+|
|  |  |   Yes No   |
|--+--+|
|Yes   |Yes   |   911 538  |
|--+--+|
|  |No|   44 456   |
|--+--+|
|No|Yes   |3 43|
|--+--+|
|  |No|2 279   |
|--+--+|



Coding and Models

  table8.3 = read.table(textConnection("alc cig mar count
  Yes Yes Yes 911
  Yes Yes No  538
  Yes No  Yes 44
  Yes No  No  456
  No  Yes Yes 3
  No  Yes No  43
  No  No  Yes 2
  No  No  No  279"),header=TRUE)
  closeAllConnections()
  # independence model (A,C,M)
  fit1.a.c.m = glm(count ~ mar+cig+alc, family=poisson, data=table8.3)
  fit1.glm$fitted.values
  # intermediate model
  fit2.m.ca = glm(count ~ mar+cig:alc, family=poisson, data=table8.3)
  fit2.m.ca$fitted.values
  # homogeneous association model
  fit3.m.c.a  =  glm(count  ~  mar:cig+mar:alc+cig:alc, family=poisson,
  data=table8.3)
  fit3.m.c.a$fitted.values
  # saturated model
  fits = glm(count ~ mar*cig*alc, family=poisson, data=table8.3)
  fits$fitted.values


The  coding  for  variables  in the above program and the fitted values are

given  below – they show that the homogeneous association model is the only

model that fits these data well.

|-+++++-+---|
| Alcohol |  Cigarette |  Marijuana | Actual |   (A,C,M)  |  (AC,M) | 
(AC:AM:CM)|
|   use   | Use| Use|  (ACM) | Independenc| | 
homogeneou|
| ||||  e | | s 
|
|-+++++-+---|
|   Yes   | Yes| Yes|   911  |540.0   |  611.2  |   910.4   
|
|-+++++-+---|
| || No |   538  |740.2   |  837.8  |   538.6   
|
|-+++++-+---|
| | No | Yes|   44   |282.1   |  210.9  |44.6   
|
|-+++++-+---|
| || No |   456  |386.7   |  289.1  |   455.4   
|
|-+++++-+---|
|No   | Yes| Yes|3   |90.6|   19.4  |3.6
|
|-+++++-+---|
| || No |   43   |124.2   |   26.6  |42.4   
|
|-+++++-+---|
| | No | Yes|2   |47.3|  118.5  |1.4
|
|-+++++-+---|
| || No |   279  |64.9|  162.5  |   279.6   
|
|-+++++-+---|







   
 Koen Hufkens  
  To 
 Sent by:  r-help

[R] how to estimate multidimensional spectral measure of coherence

2009-03-23 Thread mauede
Please, does anyone know of an R packge to estimate multidimensional spectral 
measure of coherence within a moving time window ?
Some time ago I expeimented with a similar package that performs Cross Spectrum 
Analysis on the whole signal though.
Unluckily I deal with non-stationary signals whose properties change along with 
time. Therefore estimates can only be made over time periods 
roughly proportional to the reciprocal of the rate at which properties are 
changing.
Thank you very much.
Maura


tutti i telefonini TIM!


[[alternative HTML version deleted]]

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


Re: [R] Problems with combining plots

2009-03-23 Thread baptiste auguie


On 23 Mar 2009, at 11:52, johnhj wrote:



I have still the same problem... As you said I tried with  
par(mfrow=c(2,1))
and par(mfrow=c(1,2)) but without success. Could R compiler be the  
problem ?




Why not? But may I suggest you try first the following example I sent  
you yesterday,




png()

par(mfrow=c(2,1) )

plot(1:10)

plot(1:10)

dev.off()



Does this not produce two plots on a single page?


baptiste


PS: please do read the posting guide


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






Wills, Kellie wrote:


par(mfrow=c(1,1)) will give you just one panel.  Try  
par(mfrow=c(2,1)) or

par(mfrow=c(1,2)).


-Original Message-
From: r-help-boun...@r-project.org on behalf of johnhj
Sent: Sun 3/22/2009 10:50 AM
To: r-help@r-project.org
Subject: [R]  Problems with combining plots


Hii,

I will combine some plots. Like this example here
http://www.statmethods.net/advgraphs/layout.html I tired to do it  
for 2

plots but without success.

Here is my code:

test<-read.table(file="D:/file.txt")
space<-read.table(file="D:/space.txt")

space$gruppe <- 502*rep(1:6, each=7)
x<- c(test$V1)
y<- c(test$V2)

par(mfrow=c(1,1))

png(filename = "D:/example.png", width = 640, height = 480,  
pointsize =

12,
bg = "white",  res = NA)

boxplot(V2 ~ gruppe , data = space , col = "lightgray",boxwex=0.2)
plot(panel.first=grid(ny=NULL,nx=NULL),x,y, xlab = "Zeit(sec)", ylab
="Datenrate(MBit(sec))",ylim=c(0,40), col ="purple", type ="l", main
="combined plots",lwd=2)

dev.off()


What is the mistake in my code ?


--
View this message in context:
http://www.nabble.com/Problems-with-combining-plots-tp22646692p22646692.html
Sent from the R help mailing list archive at Nabble.com.

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

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




--
View this message in context: 
http://www.nabble.com/Problems-with-combining-plots-tp22646692p22658673.html
Sent from the R help mailing list archive at Nabble.com.

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


_

Baptiste Auguié

School of Physics
University of Exeter
Stocker Road,
Exeter, Devon,
EX4 4QL, UK

Phone: +44 1392 264187

http://newton.ex.ac.uk/research/emag

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


Re: [R] Iterative Proportional Fitting, use

2009-03-23 Thread Ravi Varadhan
Hi,

Is this what you want?

#
mat <- matrix(c(65,4,22,24,6,81,5,8,0,11,85,19,4,7,3,90),4,4)

rowmarg <- rep(100, nrow(mat))  # the row margin totals that you want

colmarg <- c(90, 120, 80, 110)  # the column margin totals that you want

newmat <- loglin( outer(rowmarg, colmarg) / sum(rowmarg), margin=list(1,2),
start=mat, fit=TRUE, eps=1.e-05, iter=100)$fit 

newmat

apply(newmat, 1, sum)

apply(newmat, 2, sum)
 

Ravi.


---

Ravi Varadhan, Ph.D.

Assistant Professor, The Center on Aging and Health

Division of Geriatric Medicine and Gerontology 

Johns Hopkins University

Ph: (410) 502-2619

Fax: (410) 614-9625

Email: rvarad...@jhmi.edu

Webpage:  http://www.jhsph.edu/agingandhealth/People/Faculty/Varadhan.html







-Original Message-
From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On
Behalf Of Koen Hufkens
Sent: Monday, March 23, 2009 9:33 AM
To: Gerard M. Keogh
Cc: r-help; r-help-boun...@r-project.org
Subject: Re: [R] Iterative Proportional Fitting, use

The data I used was just an example to work upon.

My real dataset is a confusion matrix of 24x24 (and 17x17), so coding it
into a model with all different kinds of combinations seems tedious.
That's why I hoped to use the ipf() function as it accepts a matrix as
input.

Thanks for the suggestion but I hope I can get around this without recoding
the original data.

Kind regards,
Koen

-Original Message-
From: Gerard M. Keogh 
To: Koen Hufkens 
Cc: r-help , r-help-boun...@r-project.org
Subject: Re: [R] Iterative Proportional Fitting, use
Date: Mon, 23 Mar 2009 13:11:15 +

Keon,

why not fit a loglinear independence model which as far as I know is the
same.

Gerard

Here's an example from Agresti - Intro to Cat Data analysis
Example: Alcohol, cigarette, marijuana use
|--+--+|
|  Alcohol | Cigarette|Marijuana Use   |
|  |  ||
|use   |use   ||
|--+--+|
|  |  |   Yes No   |
|--+--+|
|Yes   |Yes   |   911 538  |
|--+--+|
|  |No|   44 456   |
|--+--+|
|No|Yes   |3 43|
|--+--+|
|  |No|2 279   |
|--+--+|



Coding and Models

  table8.3 = read.table(textConnection("alc cig mar count
  Yes Yes Yes 911
  Yes Yes No  538
  Yes No  Yes 44
  Yes No  No  456
  No  Yes Yes 3
  No  Yes No  43
  No  No  Yes 2
  No  No  No  279"),header=TRUE)
  closeAllConnections()
  # independence model (A,C,M)
  fit1.a.c.m = glm(count ~ mar+cig+alc, family=poisson, data=table8.3)
  fit1.glm$fitted.values
  # intermediate model
  fit2.m.ca = glm(count ~ mar+cig:alc, family=poisson, data=table8.3)
  fit2.m.ca$fitted.values
  # homogeneous association model
  fit3.m.c.a  =  glm(count  ~  mar:cig+mar:alc+cig:alc, family=poisson,
  data=table8.3)
  fit3.m.c.a$fitted.values
  # saturated model
  fits = glm(count ~ mar*cig*alc, family=poisson, data=table8.3)
  fits$fitted.values


The  coding  for  variables  in the above program and the fitted values are

given  below – they show that the homogeneous association model is the only

model that fits these data well.

|-+++++-+---
|
| Alcohol |  Cigarette |  Marijuana | Actual |   (A,C,M)  |  (AC,M) |
(AC:AM:CM)|
|   use   | Use| Use|  (ACM) | Independenc| |
homogeneou|
| ||||  e | | s
|
|-+++++-+---
|
|   Yes   | Yes| Yes|   911  |540.0   |  611.2  |
910.4   |
|-+++++-+---
|
| || No |   538  |740.2   |  837.8  |
538.6   |
|-+++++-+---
|
| | No | Yes|   44   |282.1   |  210.9  |
44.6   |
|-++---

[R] how to save a plot in a given size in inches or centimeters

2009-03-23 Thread Lo_Lo

Hi there !

I'm ploting graphics and I'd like to save them as a .jpeg file for example,
but with a given size (in inches or cm). 

I tryed the function windows() but I think it just changes the size of the
window and not the size of the graph that you're saving.

Then I tryed with the function par(din=(width=... height=...) ) but I guess
it's protected and I can't change the values.

Finally I tryed with jpeg(file, units="in", width=..., height=...) but R
asks me for another argument "res" that is supposed to be res=72 dpi or
res=NA (I found those values in the R-help). It still doesn't work...
Because of problems of margin, or "can't create a bitmap file" or something
different again...

Could you give me some advice guys ? It sounds basic to me :) but I'm
stupidly trapped and I can't find any solutions !! 

Thanks a lot 
-- 
View this message in context: 
http://www.nabble.com/how-to-save-a-plot-in-a-given-size-in-inches-or-centimeters-tp22656478p22656478.html
Sent from the R help mailing list archive at Nabble.com.

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


Re: [R] using wavelet transform to calculate mean frequency of a signal

2009-03-23 Thread stephen sefick
global wavelet spectrum?  There is something somewhere - I just can't
remember where off the top of my head.

Stephen Sefick

On Sun, Mar 22, 2009 at 2:21 PM, stvienna wiener  wrote:
> Dear list,
>
> in short: I would like to calculate the mean frequency
> of a signal (e.g. time series) using the wavelet transform.
>
>
>
> with details: I did not find a R function to calculate a
> mean frequency using one of the cran packages.
>
> My searches using  R Site Search returned:
> **(1)waveclock {waveclock}    R Documentation
> Reconstruction of the modal frequencies in a time series using continuous
> wavelet transformation and the "crazy climbers" algorithm
> (this look like what I wanted to do, but I can't figure it out)
>
> **(2) tfmean {Rwave}    R Documentation
> Average frequency by frequency
> (this seems not to help either)
>
>
> I am a bit desperate at the moment
> and would very, very much appreciate any help.
>
> Regards,
> Stephan
>
>        [[alternative HTML version deleted]]
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>



-- 
Stephen Sefick

Let's not spend our time and resources thinking about things that are
so little or so large that all they really do for us is puff us up and
make us feel like gods.  We are mammals, and have not exhausted the
annoying little problems of being mammals.

-K. Mullis

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


Re: [R] how to save a plot in a given size in inches or centimeters

2009-03-23 Thread stephen sefick
isn't there a width height argument in the jpeg function?
?jpeg
I am probably wrong,

Stephen Sefick

On Mon, Mar 23, 2009 at 10:38 AM, Lo_Lo  wrote:
>
> Hi there !
>
> I'm ploting graphics and I'd like to save them as a .jpeg file for example,
> but with a given size (in inches or cm).
>
> I tryed the function windows() but I think it just changes the size of the
> window and not the size of the graph that you're saving.
>
> Then I tryed with the function par(din=(width=... height=...) ) but I guess
> it's protected and I can't change the values.
>
> Finally I tryed with jpeg(file, units="in", width=..., height=...) but R
> asks me for another argument "res" that is supposed to be res=72 dpi or
> res=NA (I found those values in the R-help). It still doesn't work...
> Because of problems of margin, or "can't create a bitmap file" or something
> different again...
>
> Could you give me some advice guys ? It sounds basic to me :) but I'm
> stupidly trapped and I can't find any solutions !!
>
> Thanks a lot
> --
> View this message in context: 
> http://www.nabble.com/how-to-save-a-plot-in-a-given-size-in-inches-or-centimeters-tp22656478p22656478.html
> Sent from the R help mailing list archive at Nabble.com.
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>



-- 
Stephen Sefick

Let's not spend our time and resources thinking about things that are
so little or so large that all they really do for us is puff us up and
make us feel like gods.  We are mammals, and have not exhausted the
annoying little problems of being mammals.

-K. Mullis

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


[R] Problems with adapt

2009-03-23 Thread Cristian Angelo Guevara
Hi:

I'm trying to estimate a model which involves the estimation of double
integrals, so I'm using adapt procedure. I need to calculate the
integrals trough my 2000 size database, so I do it using a loop. My
code correctly estimates the integral for the first row, but for the
second R crashes. I tried changing the order of the data but the
result is the same so I guess there is something wrong with my code of
the loop or, more probably, on the use of adapt. I'm new to R so I
would appreciate any comment on my code below.

Thanks
Angelo


--

DATA<-read.table("Data.csv",header=TRUE,sep=",")
library(adapt)

mnl.lik<-function(theta,y){
th1<-theta[1]
th2<-theta[2]
tha<-theta[3]
thb<-theta[4]
thc<-theta[5]
thp<-theta[6]
thmu<-theta[7]
alfz1<-theta[8]
alfz2<-theta[9]
alfc<-theta[10]
alf<-theta[11]
r<- 1
s<- 1
n<-2000
lik<-numeric(n)
int<-numeric(n)
v<- numeric(2)

for (i in 1:n) {  #Beggin Loop

lstarpre<- function(v){   #This is the fuction to be integrated

e1<-y$p1[i]-alfz1*y$z_a1[i] - alfz2*y$z_b1[i] -alfc*y$c1[i] -alf
e2<-y$p2[i]-alfz1*y$z_a2[i] - alfz2*y$z_b2[i] -alfc*y$c2[i] -alf
e3<-y$p3[i]-alfz1*y$z_a3[i] - alfz2*y$z_b3[i] -alfc*y$c3[i] -alf

U1<- (th1 +tha*y$a1[i] +thb*y$b1[i] +thc*y$c1[i]  +thp*y$p1[i]
+thmu*(e1[i]+v[1]))
U2<- (th2 +tha*y$a2[i] +thb*y$b2[i] +thc*y$c2[i]  +thp*y$p2[i]
+thmu*(e2[i]+v[2]))
U3<- (+tha*y$a3[i] +thb*y$b3[i] +thc*y$c3[i]  +thp*y$p3[i] +thmu*(e3[i]))

Pni<- (y$Ch1[i]*exp(U1) + y$Ch2[i]*exp(U2)+
y$Ch3[i]*exp(U3))/(exp(U1)+exp(U2)+exp(U3))
fe<- exp((-e1[i]*e1[i]-e2[i]*e2[i]-e3[i]*e3[i])/(2*r))

fv<- exp((-v[1]*v[1]-v[2]*v[2])/(2*s))
prelik<- Pni*fe*fv
return(prelik)
}

int[i]<- adapt(2,lo = c(-4,-4), up = c(4,4), functn= lstarpre, eps=0.01)$value
print(int[i])
lik[i]<- log(max(1E-20,int))
print(lik[i])
}   #End Loop

logl<-sum(lik)
return(-logl)
}

p<-optim(c(0,0,0,0,0,0,0,0,0,0,0),mnl.lik,y=DATA,method="BFGS",hessian=T,
control = list(maxit=6000, temp=2000, trace=4))

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


Re: [R] how to estimate multidimensional spectral measure of coherence

2009-03-23 Thread stephen sefick
have you tried sowas?  I know you had talked about it, but it may do
what you want.  I have used it for the wavelet cross spectrum.

On Mon, Mar 23, 2009 at 9:47 AM,   wrote:
> Please, does anyone know of an R packge to estimate multidimensional spectral 
> measure of coherence within a moving time window ?
> Some time ago I expeimented with a similar package that performs Cross 
> Spectrum Analysis on the whole signal though.
> Unluckily I deal with non-stationary signals whose properties change along 
> with time. Therefore estimates can only be made over time periods
> roughly proportional to the reciprocal of the rate at which properties are 
> changing.
> Thank you very much.
> Maura
>
>
> tutti i telefonini TIM!
>
>
>        [[alternative HTML version deleted]]
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>



-- 
Stephen Sefick

Let's not spend our time and resources thinking about things that are
so little or so large that all they really do for us is puff us up and
make us feel like gods.  We are mammals, and have not exhausted the
annoying little problems of being mammals.

-K. Mullis

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


Re: [R] 'require' equivalent for local functions

2009-03-23 Thread JiHO
Thanks very much to everyone. I think I will use a combination of both  
techniques.


On 2009-March-22  , at 20:08 , Duncan Murdoch wrote:

That's pretty hard to make bulletproof.  Why not just put those  
functions in a package, and use that package?



I know it will be impossible to make bullet proof and efficient at the  
same time. However, my functions are pretty specific to each project,  
have long names and do not collide with variable names (because I use  
dots in function names but camel case in variable names) so just  
looking for the name should be OK. Plus I have a simple keyboard  
shortcut in my text editor to source the current file in the currently  
R session, so it will be easy to re-source some files after I modify  
them.


On the other hand I have a bundle of general enough functions that I  
import in many projects (http://github.com/jiho/r-utils/ for those  
that might be interested in netCDF data handling, 2D arrow fields and  
ggplot2 stuff) and this one is a good candidate to be turned into a  
package.


So thanks again to everyone. Sincerely,

JiHO
---
http://jo.irisson.free.fr/

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


[R] Scaled MPSE as a test for regressors?

2009-03-23 Thread Zhou Fang
Hi,

This is really more a stats question than a R one, but

Does anyone have any familiarity with using the mean prediction
squared error scaled by the variance of the response, as a 'scale
free' criterion for evaluating different regression algorithms.

E.g.

Generate X_train, Y_train, X_test, Y_test from true f. X_test/Y_test
are generated without noise, maybe?

Use X_train, Y_train and the algorithm to make \hat{f}

Look at var(Y_test - \hat{f}(X_test))/var(Y_test)

(Some of these var maybe should be replaced with mean squared values instead.)


It seems sort of reasonable to me. You get a number between zero and
one out of it, with 1 the solution for constant fits. Anyone seen
anything like this, or know anything about properties? Has it got a
name?

Zhou Fang

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


[R] R: how to estimate multidimensional spectral measure of coherence

2009-03-23 Thread mauede
Does it work on a sliding window ? Does it estimate the cospectrum (the real 
part) and the quadrature spectrum (complex), the coherence squared, and the 
phase difference between two vector time series ?
SOme time ago I started to read its author's thesis. It seemed to me strictly 
tailored on the specific problems he studied more than a wider purpose toolkit.

Maura 

have you tried sowas?  I know you had talked about it, but it may do
what you want.  I have used it for the wavelet cross spectrum.

On Mon, Mar 23, 2009 at 9:47 AM,   wrote:
> Please, does anyone know of an R packge to estimate multidimensional spectral 
> measure of coherence within a moving time window ?
> Some time ago I expeimented with a similar package that performs Cross 
> Spectrum Analysis on the whole signal though.
> Unluckily I deal with non-stationary signals whose properties change along 
> with time. Therefore estimates can only be made over time periods
> roughly proportional to the reciprocal of the rate at which properties are 
> changing.
> Thank you very much.
> Maura
>
>
> tutti i telefonini TIM!
>
>
>        [[alternative HTML version deleted]]
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>



-- 
Stephen Sefick

Let's not spend our time and resources thinking about things that are
so little or so large that all they really do for us is puff us up and
make us feel like gods.  We are mammals, and have not exhausted the
annoying little problems of being mammals.

-K. Mullis




tutti i telefonini TIM!


[[alternative HTML version deleted]]

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


[R] aggregate() example fails

2009-03-23 Thread Kenneth Roy Cabrera Torres
Hi R users and developers.

I compile the 
R version 2.8.1 Patched (2009-03-18 r48193)
On my UBUNTU linux distribution.

But hen I ask for the aggregate example it fails.

What am I missing?

example(aggregate)

aggrgt> ## Compute the averages for the variables in 'state.x77',
grouped
aggrgt> ## according to the region (Northeast, South, North Central,
West) that
aggrgt> ## each state belongs to.
aggrgt> aggregate(state.x77, list(Region = state.region), mean)
Error en FUN(X[[1L]], ...) : 
  elemento 1 esta vacio;
  la parte de la lista de argumentos  'is.list' ha sido evaluada:
   (INDEX)

I discovery this behavior when I try to use old scripts that uses
the aggregate() function, but now it fails on my old scripts.

Thank you for your help.

Kenneth 
PD: R.Version() returns:
$platform
[1] "x86_64-unknown-linux-gnu"
$arch
[1] "x86_64"
$os
[1] "linux-gnu"
$syste
[1] "x86_64, linux-gnu"
$status
[1] "Patched"
$major
[1] "2"
$minor
[1] "8.1"
$year
[1] "2009"
$month
[1] "03"
$day
[1] "18"
$`svn rev`
[1] "48193"
$language
[1] "R"
$version.string
[1] "R version 2.8.1 Patched (2009-03-18 r48193)"

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


Re: [R] R: Multi-line texts in plots

2009-03-23 Thread Greg Snow
Look at the textplot function in the gplots package and the addtable2plot 
function in the plotrix package.  Either of those make it easy to 'plot' a 
matrix.

Hope this helps,

-- 
Gregory (Greg) L. Snow Ph.D.
Statistical Data Center
Intermountain Healthcare
greg.s...@imail.org
801.408.8111


> -Original Message-
> From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-
> project.org] On Behalf Of mau...@alice.it
> Sent: Friday, March 20, 2009 11:41 PM
> To: baptiste auguie
> Cc: r-help@r-project.org
> Subject: [R] R: Multi-line texts in plots
> 
> Now that I have my list of flags with theri respective values (thanks
> to all those who posted their suggestions):
> 
>Flags Values
> 1   TrendOff  0
> 2  MOdwt  1
> 3ZeroPadding  1
> 4 Step1HSOff  1
> 5  Step1NumHSOff  4
> 6  Step1NumHSOff  1
> 7Step1ExtOff  1
> 8 Step2HSOff  1
> 9  Step2NumHSOff  4
> 10 Step2NumHSOff  1
> 11   Step2ExtOff  1
> 12Step3HSOff  1
> 13 Step3NumHSOff  4
> 14 Step3NumHSOff  1
> 15   Step3ExtOff  1
> 16Step4HSOff  1
> 17 Step4NumHSOff  4
> 18 Step4NumHSOff  1
> 19   Step4ExtOff  1
> 
> The next step is to insert the above list in a composite plot made up
> of 4 quadrants. The top 2 quadrants are filled with a barplot, the
> bottom left quadrant is filled with a 2-tracks plot. I would like to
> place the flags list in the bottom right quadrant.
> Hiow can I do that ? Maybe through using the text command recursively
> (once for each flag pairlist) ?
> Please, notice taht what none of the graphs represent the flags value.
> Therefore a legend is unappropriate here.
> Here is a sketched incomplete and undebugged script to generate the
> composite plot I mentioned:
> 
>  alpha <- matrix(nrow=21,ncol=2,byrow=T)
>  zm<- matrix(nrow=10,ncol=2,byrow=T)
>  colnames(alpha) <- c("MEM","SpAn")
>  colnames(zm) <- c("MEM","SpAn")
> 
>  alpha[,1] <- c( 0.9453125,
>  1,
>  1,
>  1,
>  0.9765625,
>  0.9973958,
>  0.9375,
>  0.9921671,
>  0.9765625,
>  0.9322917,
>  0.96875,
>  0.8645833,
>  0.8723958,
>  0.9270833,
>  0.9791667,
>  0.9869792,
>  0.9661458,
>  0.9895833,
>  1,
>  1,
>  0.9713542 )
> 
>  alpha[,2] <- c( 0.7265625,
>  0.8828125,
>  0.7734375,
>  0.7734375,
>  0.6875,
>  0.8359375,
>  0.6484375,
>  0.8046875,
>  0.71875,
>  0.734375,
>  0.71875,
>  0.7265625,
>  0.6796875,
>  0.7265625,
>  0.6953125,
>  0.6953125,
>  0.703125,
>  0.7890625,
>  0.828125,
>  0.8359375,
>  0.8515625 )
> 
>  zm[,1] <- c( 0.857, 0, 0, 0, 0.0476, 0, 0.0476, 0, 0, 0.0476 )
> 
>  zm[,2] <- c( 0, 0, 0, 0.762, 0, 0.0476, 0.1904, 0, 0, 0 )
> 
>  flags_val <- c( 0,1,1,1,4,1,1,1,4,1,1,1,4,1,1,1,4,1,1)
>  flags_nam <- c("TrendOff","MOdwt","ZeroPadding",
> 
> "Step1HSOff","Step1NumHSOff","Step1NumHSOff","Step1ExtOff",
> 
> "Step2HSOff","Step2NumHSOff","Step2NumHSOff","Step2ExtOff",
> 
> "Step3HSOff","Step3NumHSOff","Step3NumHSOff","Step3ExtOff",
> 
> "Step4HSOff","Step4NumHSOff","Step4NumHSOff","Step4ExtOff" )
> 
>  flags <- data.frame(Flags=flags_nam,
> Values=flags_val,stringsAsFactors=FALSE)
> 
>  par(mfrow=c(2,2))
> 
>  barplot(zm[,1],width=1,space=0,horiz=FALSE,axes=FALSE,main="MEM
> Vanishing Moments Distribution",cex.main=0.95)
>  xxLab <- seq(0,11,1)
>  axis(1,at=xxLab-
> 0.5,labels=as.character(xxLab),col="red",col.axis="red",font.axis=1,xpd
> =TRUE,
>   cex.lab=1)
>  yyLab <- seq(0,1,0.14)
> 
> axis(2,at=yyLab,labels=as.character(yyLab),col="red",col.axis="red",fon
> t.axis=1,xpd=TRUE,
>   cex.lab=0.8,cex.axis=0.8)
>  barplot(zm[,2],width=1,space=0,horiz=TRUE,axes=FALSE)
>  text(x=25.5,y=15.6,pos=4,"SpAn Vanishing Moments
> Distribution",cex=1.2,font=2)
>  axis(2,at=xxLab-
> 1,labels=as.character(xxLab),col="red",col.axis="red",font.axis=1,xpd=T
> RUE,
>   cex.lab=1)
>  axis(3,at=yyLab-
> 1,labels=as.character(yyLab),col="red",col.axis="red",font.axis=1,xpd=T
> RUE,
>   cex.lab=0.8,cex.axis=0.8)
>  cycle <- seq(1,21,1)
>  plot(cycle, alpha[,1], pch=1, type = "l", col=3, lty=2, main="Alpha")
>  points(cycle, alpha[,2], pch=3, type="b", lty=1, pch=4, col=6)
>  legend("topright", legend=c("MEM","SpAn"),
> col=c(3,6),text.col="green4",lty=c(2,1), pch=c(1,3),merge=FALSE,
> bg='gray90')
> 
> Thank tyou so much.
> Maura
> 
> 
> 
> 
> 
> 
> 
> 
> 
> 
> -Messaggio originale-
> Da: bapti

Re: [R] how to save a plot in a given size in inches or centimeters

2009-03-23 Thread David M Smith
On Mon, Mar 23, 2009 at 7:38 AM, Lo_Lo  wrote:
> I'm ploting graphics and I'd like to save them as a .jpeg file for example,
> but with a given size (in inches or cm).
>
> I tryed the function windows() but I think it just changes the size of the
> window and not the size of the graph that you're saving.
>
> Then I tryed with the function par(din=(width=... height=...) ) but I guess
> it's protected and I can't change the values.
>
> Finally I tryed with jpeg(file, units="in", width=..., height=...) but R
> asks me for another argument "res" that is supposed to be res=72 dpi or
> res=NA (I found those values in the R-help). It still doesn't work...
> Because of problems of margin, or "can't create a bitmap file" or something
> different again...
>
> Could you give me some advice guys ? It sounds basic to me :) but I'm
> stupidly trapped and I can't find any solutions !!
>
> Thanks a lot

You'll be wanting to call the jpeg() function directly with the width=
and height= arguments.  These are given in pixels, which makes things
a little tricky if you want to measure in inches -- you'll need also
use the "res" argument to set the resolution.  See:

http://blog.revolution-computing.com/2009/01/10-tips-for-making-your-r-graphics-look-their-best.html

for some examples.  I'm guessing you want to print your graphic in
which case jpeg() probably isn't the best choice (it rarely is for
statistical graphics).  See that same post for some other suggestions.

# David Smith

-- 
David M Smith 
Director of Community, REvolution Computing www.revolution-computing.com
Tel: +1 (206) 577-4778 x3203 (Seattle, USA)

Check out our upcoming events schedule at www.revolution-computing.com/events

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


Re: [R] variance/mean

2009-03-23 Thread Bert Gunter
Inline Below.

-- Bert 


Bert Gunter
Genentech Nonclinical Biostatistics
650-467-7374

-Original Message-
From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On
Behalf Of Wacek Kusnierczyk
Sent: Sunday, March 22, 2009 2:16 AM
To: rkevinbur...@charter.net
Cc: r-help@r-project.org
Subject: Re: [R] variance/mean

rkevinbur...@charter.net wrote:
> At the risk of appearing ignorant why is the folowing true?
>
> o <- cbind(rep(1,3),rep(2,3),rep(3,3))
> var(o)
>  [,1] [,2] [,3]
> [1,]000
> [2,]000
> [3,]000
>
> and
>
> mean(o)
> [1] 2
>
> How do I get mean to return an array similar to var? I would expect in the
above example a vector of length 3 {1,2,3}.
>   
You said:

"you may well be ignorant about how var works with matrices, but this
does not mean it's your fault.  the documentation is typically cryptical."


-- How so? ?var clearly states:

" ... If x and y are matrices then the covariances (or correlations) between
the columns of x and the columns of y are computed. "

and the Arguments section says:

x a numeric vector, matrix or data frame. 
y NULL (default) or a vector, matrix or data frame with compatible
dimensions to x. The default is equivalent to y = x (but more efficient). 


This is as clear as I would know how to state. I think "...typically
cryptical" is a canard and most unfair.

-- Bert

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


[R] Looping of read.table and assignment

2009-03-23 Thread Steve Murray

Dear all,

I am trying to read in and assign data from 50 tables in an automated fashion. 
I have the following code, which I created with the help of textbooks and the 
internet, but it only seems to read in the final data file over and over again. 
For example, when I type:> table_1951  I get the same values in the table as 
when I type> table_2000  despite the values in the source tables being 
different:


year <- 1951:2000
filelist <- paste("C:\\Documents and Settings\\Data\\table_",year,".txt", 
sep="")
filelist
# Code seems to operate successfully up to this point

for (i in filelist) {
   for (iyear in 1951:2000) {
  assign(paste("table_",iyear, sep=""),read.table(file=i, header=TRUE, 
sep=","))
  noquote(paste("LOADED FILE:",paste("table_",iyear, sep=""),sep=" "))
 }
 }


Can anyone see what I've done wrong here?

And just as an aside, as you can see, I've inserted the 'noquote' line so that 
when the code is running I should be able to see each file being read in - 
mainly as a 'checker'. Should this work as anticipated, with each line being 
displayed with its corresponding table number after it's been read in?

Many thanks for any help offered,

Steve


_
[[elided Hotmail spam]]

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


Re: [R] Plot and Boxplot in the same graph

2009-03-23 Thread Greg Snow
Look at the symbols function or the subplot function in the TeachingDemos 
package (or my.symbols in TeachingDemos).

Hope this helps,

-- 
Gregory (Greg) L. Snow Ph.D.
Statistical Data Center
Intermountain Healthcare
greg.s...@imail.org
801.408.8111

> -Original Message-
> From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-
> project.org] On Behalf Of johnhj
> Sent: Sunday, March 22, 2009 5:18 AM
> To: r-help@r-project.org
> Subject: Re: [R] Plot and Boxplot in the same graph
> 
> 
> I tried to do it with:
> 
> > y <- rnorm(100)
> > x <- gl(2,50)
> > boxplot(x,y)
> > points(x,y)
> 
> But the problem is, that the the y coordinates are shown for the
> boxplot and
> not for "points(x,y)"
> Is it possible to show the graph with the (x,y) coordinates with the
> points() function and the boxplots only for the x coordinates. A better
> solution could be to have a seperated y axis on the left side for the
> boxplots().
> 
> Is it possible to do it in this was ?
> 
> 
> 
> 
> Paul Johnson-11 wrote:
> >
> > On Fri, Mar 20, 2009 at 10:02 PM, johnhj  wrote:
> >>
> >> Hii,
> >>
> >> Is it possible, to use the plot() funktion and the boxplot()
> funktion
> >> together ?
> >> I will plot a simple graph and additionally to the graph on certain
> >> places
> >> boxplots. I have imagined to plot the graph a little bit
> transparency and
> >> show in the same graph on certain places boxplots
> >>
> >> Is it possible to do it in this way ?
> >>
> >> greetings,
> >> johnh
> >> --
> >>
> > Run the boxplot first, then use points() or other subsidiary plot
> > functions to add the points in the figure.
> >
> >> y <- rnorm(100)
> >> x <- gl(2,50)
> >> boxplot(x,y)
> >> points(x,y)
> >
> >
> > --
> > Paul E. Johnson
> > Professor, Political Science
> > 1541 Lilac Lane, Room 504
> > University of Kansas
> >
> > __
> > R-help@r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide
> > http://www.R-project.org/posting-guide.html
> > and provide commented, minimal, self-contained, reproducible code.
> >
> >
> 
> --
> View this message in context: http://www.nabble.com/Plot-and-Boxplot-
> in-the-same-graph-tp22632355p22645076.html
> Sent from the R help mailing list archive at Nabble.com.
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-
> guide.html
> and provide commented, minimal, self-contained, reproducible code.

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


Re: [R] Looping of read.table and assignment

2009-03-23 Thread Gabor Grothendieck
See:

https://stat.ethz.ch/pipermail/r-help/2009-March/192422.html

On Mon, Mar 23, 2009 at 12:52 PM, Steve Murray  wrote:
>
> Dear all,
>
> I am trying to read in and assign data from 50 tables in an automated 
> fashion. I have the following code, which I created with the help of 
> textbooks and the internet, but it only seems to read in the final data file 
> over and over again. For example, when I type:> table_1951  I get the same 
> values in the table as when I type> table_2000  despite the values in the 
> source tables being different:
>
>
> year <- 1951:2000
> filelist <- paste("C:\\Documents and Settings\\Data\\table_",year,".txt", 
> sep="")
> filelist
> # Code seems to operate successfully up to this point
>
> for (i in filelist) {
>   for (iyear in 1951:2000) {
>          assign(paste("table_",iyear, sep=""),read.table(file=i, header=TRUE, 
> sep=","))
>          noquote(paste("LOADED FILE:",paste("table_",iyear, sep=""),sep=" "))
>                     }
>                     }
>
>
> Can anyone see what I've done wrong here?
>
> And just as an aside, as you can see, I've inserted the 'noquote' line so 
> that when the code is running I should be able to see each file being read in 
> - mainly as a 'checker'. Should this work as anticipated, with each line 
> being displayed with its corresponding table number after it's been read in?
>
> Many thanks for any help offered,
>
> Steve
>
>
> _
> [[elided Hotmail spam]]
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

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


[R] Basic regression output question

2009-03-23 Thread GRANT Lewis
Hi

Probably a very basic question: 

I am regressing a matrix of 50 response variables against a matrix of 10
factors using the lm function. This gives me an object with the output
for 50 regressions, as required. How do I now "use" the data? 

For example, below is some code that generates sample data and shows me
the summary statistics for regression 50. If I want to access the data,
perhaps to extract a particular column or row, how do I do this? Code
such as (coef(summary(lm(t(returns)~factors)))[50])[1,4], which works
when I have only one response variable, returns me an error.

Thanks in advance 

Lewis


EXAMPLE CODE

factors<-matrix(runif(400),nrow=40)
returns1<-matrix(runif(40),nrow=1)
returns2<-matrix(runif(2000),nrow=50)
coef(summary(lm(t(returns1)~factors)))[1,4]

[1] 0.01062590

(coef(summary(lm(t(returns2)~factors)))[50])

Response Y50 :
   Estimate Std. Error t value  Pr(>|t|)
(Intercept)  0.54829326  0.3230444  1.69726923 0.1003545
factors1 0.08225028  0.1858496  0.44256369 0.6613648
factors2-0.03084209  0.1968324 -0.15669216 0.8765733
factors3-0.11287874  0.1845361 -0.61168926 0.5455093
factors4 0.01153631  0.2197322  0.05250168 0.9584889
factors5-0.21893452  0.2079244 -1.05295235 0.3010575
factors6 0.12978584  0.2164691  0.59955814 0.5534567
factors7 0.09461358  0.1918524  0.49315832 0.6256146
factors8 0.01592064  0.2096206  0.07594975 0.9399806
factors9-0.06970911  0.2072088 -0.33641968 0.7389767
factors10   -0.04185789  0.1910195 -0.21912894 0.8280846

(coef(summary(lm(t(returns2)~factors)))[50])[1,4]

Error in `[.default`((coef(summary(lm(t(returns2) ~ factors)))[50]), 1,
: 
  incorrect number of dimensions





**
Hermes Fund Managers Limited 
Registered in England No. 1661776, Lloyds Chambers, 1 Portsoken Street, London 
E1 8HZ

*** Please read the Hermes email disclaimer at 
http://www.hermes.co.uk/email_terms.htm before acting on this email or opening 
any attachment ***

The contents of this email are confidential.  If you have received this message 
in error, please delete it immediately and contact the sender directly or the 
Hermes IT Helpdesk on +44(0)20 7680 2117.  Any reliance on, use, disclosure, 
dissemination, distribution or copying of this email is unauthorised and 
strictly prohibited.

This message has been checked for viruses but the recipient is strongly advised 
to rescan the message before opening any attachments or attached executable 
files.  Hermes do not accept any liability for any damage sustained as a result 
of a virus introduced by this email or any attachment.


**


__
This email has been scanned by the MessageLabs Email Security System.
For more information please visit http://www.messagelabs.com/email 

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


Re: [R] Looping of read.table and assignment

2009-03-23 Thread baptiste auguie

You need only one loop,


year <- 1951:2000
filelist <- paste("C:\\Documents and Settings\\Data\ 
\table_",year,".txt", sep="")

filelist



for (ii in seq_along(year)) {
 assign(paste("table_", year[ii], sep=""),
read.table(file=ifile[ii], header=TRUE, 
sep=","))
}




or better yet,


d <- lapply( filelist, read.csv)
names(d) <- paste(year)



head(d["1952"])



(untested)

HTH,

baptiste

On 23 Mar 2009, at 16:52, Steve Murray wrote:



Dear all,

I am trying to read in and assign data from 50 tables in an  
automated fashion. I have the following code, which I created with  
the help of textbooks and the internet, but it only seems to read in  
the final data file over and over again. For example, when I type:>  
table_1951  I get the same values in the table as when I type>  
table_2000  despite the values in the source tables being different:



year <- 1951:2000
filelist <- paste("C:\\Documents and Settings\\Data\ 
\table_",year,".txt", sep="")

filelist
# Code seems to operate successfully up to this point

for (i in filelist) {
  for (iyear in 1951:2000) {
 assign(paste("table_",iyear, sep=""),read.table(file=i,  
header=TRUE, sep=","))
 noquote(paste("LOADED FILE:",paste("table_",iyear,  
sep=""),sep=" "))

}
}


Can anyone see what I've done wrong here?

And just as an aside, as you can see, I've inserted the 'noquote'  
line so that when the code is running I should be able to see each  
file being read in - mainly as a 'checker'. Should this work as  
anticipated, with each line being displayed with its corresponding  
table number after it's been read in?


Many thanks for any help offered,

Steve


_
[[elided Hotmail spam]]

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


_

Baptiste Auguié

School of Physics
University of Exeter
Stocker Road,
Exeter, Devon,
EX4 4QL, UK

Phone: +44 1392 264187

http://newton.ex.ac.uk/research/emag

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


[R] Fitting multiple Gaussians to data

2009-03-23 Thread Blake Sweeney
Hello,

I am trying to fit several Gaussian distributions to a data set. I
need to know the mean, standard deviation and maximum of each
Gaussian. My data set is simply a list of signal strengths and no
other information. I do not know the number of Guassians to fit. I've
looked at some packages like mda, mfp, and sm, but I don't think these
do what I want.

Does anyone have any suggestions?


Thank you,
Blake

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


[R] Extracting SD of random effects from lme object

2009-03-23 Thread Ben Domingue
Hello,
How do I get the standard deviations for the random effects out of the
lme object?  I feel like there's probably a simple way of doing this,
but I can't see it.  Using the first example from the documentation:

> fm1 <- lme(distance ~ age, data = Orthodont) # random is ~ age
> fm1
Linear mixed-effects model fit by REML
  Data: Orthodont
  Log-restricted-likelihood: -221.3183
  Fixed: distance ~ age
(Intercept) age
 16.761   0.6601852

Random effects:
 Formula: ~age | Subject
 Structure: General positive-definite
StdDevCorr
(Intercept) 2.3270339 (Intr)
age 0.2264276 -0.609
Residual1.3100399

Number of Observations: 108
Number of Groups: 27

I want to extract the column vector (2.3270339, 0.2264276,
1.3100399)'.  Any thoughts?
Thanks,
Ben

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


[R] How to set up a function for "Central Limit Theorem"

2009-03-23 Thread pfc_ivan

Hello guys, I am stuck here:

How do I make 1000 samples of n = 10 observations from an Exponential
distribution and then compute the mean for all those 1000 samples? 

Basically I need to prove the Central Limit theorem, which states:

http://www.nabble.com/file/p22664113/d175f06cbf200bd52a2c27a2e56dc594.png 

Where the Sn is sum of random variables, n we have from the question, mu is
mean and (sigma)^2 is variance. 

I am having trouble setting up the function to do this. 

Any help apreciated!

-- 
View this message in context: 
http://www.nabble.com/How-to-set-up-a-function-for-%22Central-Limit-Theorem%22-tp22664113p22664113.html
Sent from the R help mailing list archive at Nabble.com.

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


[R] Replacing a few variable values within a DataFrame...

2009-03-23 Thread Jason Rupert

I would like to replace a few varaibles within a data frame. 

For example, in the dataframe below (contrived) I would like to replace the 
current housesize value only if the Location is HSV.   However, I would like to 
leave the other values intact.  


I tried "ifelse", but I don't really need the else condition. 

test_data2_df<-data.frame(Variables=c("SQR Footage","SQR Footage","SQR 
Footage","SQR Footage","SQR Footage","SQR Footage","SQR Footage","SQR 
Footage","SQR Footage","SQR Footage","SQR Footage","SQR Footage","SQR 
Footage","SQR Footage","SQR Footage"),HouseSize=c(10, 20, 30, 40, 50, 15, 25, 
35, 45, 55, 18, 28, 38, 48, 58), Lot=c(1, 2, 3, 4, 5, 10, 20, 30, 40, 50, 11, 
21, 31, 41, 51), Location=c("HSV", "ATH", "HSV", "ATH", "FLO", 
 "HSV", "ATH", "HSV", "ATH", "FLO", 
 "HSV", "ATH", "HSV", "ATH", "FLO"))

Moreover, I want to preserve the rest of the values within the data.  I thought 
some combination using the subset would get me there, but I ended up with a 
long and ungangly looking code. 

Thanks for any help that is offered.

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


Re: [R] Replacing a few variable values within a DataFrame...

2009-03-23 Thread Jorge Ivan Velez
Dear Jason,
Try this:

index<- with(test_data2_df,Location %in% "HSV")
test_data2_df$HouseSize<-ifelse(index,1000,test_data2_df$HouseSize)  # I
used 1000 as reference
test_data2_df

HTH,

Jorge


On Mon, Mar 23, 2009 at 1:39 PM, Jason Rupert wrote:

>
> I would like to replace a few varaibles within a data frame.
>
> For example, in the dataframe below (contrived) I would like to replace the
> current housesize value only if the Location is HSV.   However, I would like
> to leave the other values intact.
>
>
> I tried "ifelse", but I don't really need the else condition.
>
> test_data2_df<-data.frame(Variables=c("SQR Footage","SQR Footage","SQR
> Footage","SQR Footage","SQR Footage","SQR Footage","SQR Footage","SQR
> Footage","SQR Footage","SQR Footage","SQR Footage","SQR Footage","SQR
> Footage","SQR Footage","SQR Footage"),HouseSize=c(10, 20, 30, 40, 50, 15,
> 25, 35, 45, 55, 18, 28, 38, 48, 58), Lot=c(1, 2, 3, 4, 5, 10, 20, 30, 40,
> 50, 11, 21, 31, 41, 51), Location=c("HSV", "ATH", "HSV", "ATH", "FLO",
> "HSV", "ATH", "HSV", "ATH", "FLO",
> "HSV", "ATH", "HSV", "ATH", "FLO"))
>
> Moreover, I want to preserve the rest of the values within the data.  I
> thought some combination using the subset would get me there, but I ended up
> with a long and ungangly looking code.
>
> Thanks for any help that is offered.
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

[[alternative HTML version deleted]]

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


Re: [R] Replacing a few variable values within a DataFrame...

2009-03-23 Thread baptiste auguie


On 23 Mar 2009, at 17:39, Jason Rupert wrote:



I would like to replace a few varaibles within a data frame.

For example, in the dataframe below (contrived) I would like to  
replace the current housesize value only if the Location is HSV.
However, I would like to leave the other values intact.




How about,

test_data2_df[test_data2_df$Location=="HSV", ] # these are the values  
to change


test_data2_df$housesize[test_data2_df$Location=="HSV"]  <- 28.3 # or  
whatever values



HTH,

baptiste



I tried "ifelse", but I don't really need the else condition.

test_data2_df<-data.frame(Variables=c("SQR Footage","SQR  
Footage","SQR Footage","SQR Footage","SQR Footage","SQR  
Footage","SQR Footage","SQR Footage","SQR Footage","SQR  
Footage","SQR Footage","SQR Footage","SQR Footage","SQR  
Footage","SQR Footage"),HouseSize=c(10, 20, 30, 40, 50, 15, 25, 35,  
45, 55, 18, 28, 38, 48, 58), Lot=c(1, 2, 3, 4, 5, 10, 20, 30, 40,  
50, 11, 21, 31, 41, 51), Location=c("HSV", "ATH", "HSV", "ATH", "FLO",

"HSV", "ATH", "HSV", "ATH", "FLO",
"HSV", "ATH", "HSV", "ATH",  
"FLO"))


Moreover, I want to preserve the rest of the values within the  
data.  I thought some combination using the subset would get me  
there, but I ended up with a long and ungangly looking code.


Thanks for any help that is offered.

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


_

Baptiste Auguié

School of Physics
University of Exeter
Stocker Road,
Exeter, Devon,
EX4 4QL, UK

Phone: +44 1392 264187

http://newton.ex.ac.uk/research/emag

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


[R] Copying a subset of a matrix

2009-03-23 Thread akintayo holder
Hi,
I have a matrix and I want to create a matrix that includes every nth column
of the original matrix. Does anyone know how I can go about doing that ? Or
if you have an idea how I can do the same thing with a data frame ? If you
can just point me towards the appropriate function I would appreciate it.
Thanks

[[alternative HTML version deleted]]

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


Re: [R] How to set up a function for "Central Limit Theorem"

2009-03-23 Thread stephen sefick
homework?

On Mon, Mar 23, 2009 at 1:30 PM, pfc_ivan  wrote:
>
> Hello guys, I am stuck here:
>
> How do I make 1000 samples of n = 10 observations from an Exponential
> distribution and then compute the mean for all those 1000 samples?
>
> Basically I need to prove the Central Limit theorem, which states:
>
> http://www.nabble.com/file/p22664113/d175f06cbf200bd52a2c27a2e56dc594.png
>
> Where the Sn is sum of random variables, n we have from the question, mu is
> mean and (sigma)^2 is variance.
>
> I am having trouble setting up the function to do this.
>
> Any help apreciated!
>
> --
> View this message in context: 
> http://www.nabble.com/How-to-set-up-a-function-for-%22Central-Limit-Theorem%22-tp22664113p22664113.html
> Sent from the R help mailing list archive at Nabble.com.
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>



-- 
Stephen Sefick

Let's not spend our time and resources thinking about things that are
so little or so large that all they really do for us is puff us up and
make us feel like gods.  We are mammals, and have not exhausted the
annoying little problems of being mammals.

-K. Mullis

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


Re: [R] How to set up a function for "Central Limit Theorem"

2009-03-23 Thread pfc_ivan

I tried using the for (i..) to make 1000 differents sets of numbers, but then
I dont know how to get the mean value of all of them... because I dont even
think 1000 different sets of numbers were made, because when i print it it
always shows me the same values, basically I dont know how to replicate the
funcion many times. 



stephen sefick wrote:
> 
> homework?
> 
> On Mon, Mar 23, 2009 at 1:30 PM, pfc_ivan  wrote:
>>
>> Hello guys, I am stuck here:
>>
>> How do I make 1000 samples of n = 10 observations from an Exponential
>> distribution and then compute the mean for all those 1000 samples?
>>
>> Basically I need to prove the Central Limit theorem, which states:
>>
>> http://www.nabble.com/file/p22664113/d175f06cbf200bd52a2c27a2e56dc594.png
>>
>> Where the Sn is sum of random variables, n we have from the question, mu
>> is
>> mean and (sigma)^2 is variance.
>>
>> I am having trouble setting up the function to do this.
>>
>> Any help apreciated!
>>
>> --
>> View this message in context:
>> http://www.nabble.com/How-to-set-up-a-function-for-%22Central-Limit-Theorem%22-tp22664113p22664113.html
>> Sent from the R help mailing list archive at Nabble.com.
>>
>> __
>> R-help@r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide
>> http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>>
> 
> 
> 
> -- 
> Stephen Sefick
> 
> Let's not spend our time and resources thinking about things that are
> so little or so large that all they really do for us is puff us up and
> make us feel like gods.  We are mammals, and have not exhausted the
> annoying little problems of being mammals.
> 
>   -K. Mullis
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 
> 

-- 
View this message in context: 
http://www.nabble.com/How-to-set-up-a-function-for-%22Central-Limit-Theorem%22-tp22664113p22666250.html
Sent from the R help mailing list archive at Nabble.com.

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


[R] repeated measures specification in lmer

2009-03-23 Thread Lawrence Hanser
Dear Colleagues,
I have what Roger Kirk (Experimental Design: Procedures for the Behavioral
Sciences, 1968) refers to as a randomized block factorial design.  The anova
table would look like this:

  df
A 3
Subj/A  103 (error term for A)
B   23
A*B69
B*Subj/A 2369 (error term for B and A*B)

Subjects are nested within A and give a response for each B.  If y is the
dependent variable, is this the correct lmer specification for the above,
where ID is the variable name for Subj:

lmer(y ~ A + B + A*B + (A|ID))

Am I barking up the right tree?  I can also fit:

aov(y ~ A + B + A*B  + ID)
then I have to do some hand calculations to use ID as the error term for A.
 The residual (really B*ID) is the correct error term for B and A*B.

Thanks,

Larry

[[alternative HTML version deleted]]

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


Re: [R] How to set up a function for "Central Limit Theorem"

2009-03-23 Thread baptiste auguie

(a <- replicate(5,rnorm(10)))
colMeans(a)

should get you started.

HTH,

baptiste


On 23 Mar 2009, at 18:29, pfc_ivan wrote:



I tried using the for (i..) to make 1000 differents sets of numbers,  
but then
I dont know how to get the mean value of all of them... because I  
dont even
think 1000 different sets of numbers were made, because when i print  
it it
always shows me the same values, basically I dont know how to  
replicate the

funcion many times.



stephen sefick wrote:


homework?

On Mon, Mar 23, 2009 at 1:30 PM, pfc_ivan   
wrote:


Hello guys, I am stuck here:

How do I make 1000 samples of n = 10 observations from an  
Exponential

distribution and then compute the mean for all those 1000 samples?

Basically I need to prove the Central Limit theorem, which states:

http://www.nabble.com/file/p22664113/d175f06cbf200bd52a2c27a2e56dc594.png

Where the Sn is sum of random variables, n we have from the  
question, mu

is
mean and (sigma)^2 is variance.

I am having trouble setting up the function to do this.

Any help apreciated!

--
View this message in context:
http://www.nabble.com/How-to-set-up-a-function-for-%22Central-Limit-Theorem%22-tp22664113p22664113.html
Sent from the R help mailing list archive at Nabble.com.

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





--
Stephen Sefick

Let's not spend our time and resources thinking about things that are
so little or so large that all they really do for us is puff us up  
and

make us feel like gods.  We are mammals, and have not exhausted the
annoying little problems of being mammals.

-K. Mullis

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




--
View this message in context: 
http://www.nabble.com/How-to-set-up-a-function-for-%22Central-Limit-Theorem%22-tp22664113p22666250.html
Sent from the R help mailing list archive at Nabble.com.

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


_

Baptiste Auguié

School of Physics
University of Exeter
Stocker Road,
Exeter, Devon,
EX4 4QL, UK

Phone: +44 1392 264187

http://newton.ex.ac.uk/research/emag

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


Re: [R] How to set up a function for "Central Limit Theorem"

2009-03-23 Thread Peter Dalgaard

pfc_ivan wrote:

I tried using the for (i..) to make 1000 differents sets of numbers, but then
I dont know how to get the mean value of all of them... because I dont even
think 1000 different sets of numbers were made, because when i print it it
always shows me the same values, basically I dont know how to replicate the
funcion many times. 


replicate() is your friend. Don't know why you get the same numbers 
multiple times though.





stephen sefick wrote:

homework?

On Mon, Mar 23, 2009 at 1:30 PM, pfc_ivan  wrote:

Hello guys, I am stuck here:

How do I make 1000 samples of n = 10 observations from an Exponential
distribution and then compute the mean for all those 1000 samples?

Basically I need to prove the Central Limit theorem, which states:

http://www.nabble.com/file/p22664113/d175f06cbf200bd52a2c27a2e56dc594.png

Where the Sn is sum of random variables, n we have from the question, mu
is
mean and (sigma)^2 is variance.

I am having trouble setting up the function to do this.

Any help apreciated!

--
View this message in context:
http://www.nabble.com/How-to-set-up-a-function-for-%22Central-Limit-Theorem%22-tp22664113p22664113.html
Sent from the R help mailing list archive at Nabble.com.

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




--
Stephen Sefick

Let's not spend our time and resources thinking about things that are
so little or so large that all they really do for us is puff us up and
make us feel like gods.  We are mammals, and have not exhausted the
annoying little problems of being mammals.

-K. Mullis

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







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

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


Re: [R] Replacing a few variable values within a DataFrame...

2009-03-23 Thread Jason Rupert

Phil and all, 

Thank you for the responses.  That method worked great!


However, I did try to replace the a few of the "Variables" using the same 
method and received an error...

test_data2_df[test_data2_df$Location == 'HSV','Variables'] = "YardSize"

Is there something special I need to do in order to replace the string 
characters? 

Right now it appears "NA" is being used. 

Thank you again.

test_data2_df<-data.frame(Variables=c("SQR Footage","SQR Footage","SQR 
Footage","SQR Footage","SQR Footage","SQR Footage","SQR Footage","SQR 
Footage","SQR Footage","SQR Footage","SQR Footage","SQR Footage","SQR 
Footage","SQR Footage","SQR Footage"), 
HouseSize=c(10, 20, 30, 40, 50, 15, 25, 35, 45, 55, 18, 28, 38, 48, 58), 
Lot=c(1, 2, 3, 4, 5, 10, 20, 30, 40, 50, 11, 21, 31, 41, 51), 
Location = c("HSV", "ATH", "HSV", "ATH", "FLO", "HSV", "ATH", "HSV", "ATH",
"FLO", "HSV", "ATH", "HSV", "ATH", "FLO"))


--- On Mon, 3/23/09, Phil Spector  wrote:

> From: Phil Spector 
> Subject: Re: [R] Replacing a few variable values within a DataFrame...
> To: "Jason Rupert" 
> Date: Monday, March 23, 2009, 12:48 PM
> Jason -
> When you say that you want to change some values (if),
> but 
> leave others the same (else) I think that the else is very 
> important.  Here's how to use ifelse to do what you
> want:
> 
> test_data2_df$HouseSize = ifelse(test_data2_df$Location
> == 'HSV',newvalue, test_data2_df$HouseSize)
> 
> Here's another way to selectively change values in a
> data frame:
> 
> test_data2_df[test_data2_df$Location ==
> 'HSV','HouseSize'] = newvalues
> 
> 
> - Phil Spector
>Statistical Computing Facility
>Department of Statistics
>UC Berkeley
>spec...@stat.berkeley.edu
> 
> 
> 
> On Mon, 23 Mar 2009, Jason Rupert wrote:
> 
> >
> > I would like to replace a few varaibles within a data
> frame.
> >
> > For example, in the dataframe below (contrived) I
> would like to replace the current housesize value only if
> the Location is HSV.   However, I would like to leave the
> other values intact.
> >
> >
> > I tried "ifelse", but I don't really
> need the else condition.
> >
> > test_data2_df<-data.frame(Variables=c("SQR
> Footage","SQR Footage","SQR
> Footage","SQR Footage","SQR
> Footage","SQR Footage","SQR
> Footage","SQR Footage","SQR
> Footage","SQR Footage","SQR
> Footage","SQR Footage","SQR
> Footage","SQR Footage","SQR
> Footage"),HouseSize=c(10, 20, 30, 40, 50, 15, 25, 35,
> 45, 55, 18, 28, 38, 48, 58), Lot=c(1, 2, 3, 4, 5, 10, 20,
> 30, 40, 50, 11, 21, 31, 41, 51), Location=c("HSV",
> "ATH", "HSV", "ATH",
> "FLO",
> > "HSV",
> "ATH", "HSV", "ATH",
> "FLO",
> > "HSV",
> "ATH", "HSV", "ATH",
> "FLO"))
> >
> > Moreover, I want to preserve the rest of the values
> within the data.  I thought some combination using the
> subset would get me there, but I ended up with a long and
> ungangly looking code.
> >
> > Thanks for any help that is offered.
> >
> > __
> > R-help@r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> > and provide commented, minimal, self-contained,
> reproducible code.
> >

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


Re: [R] Problems with combining plots

2009-03-23 Thread johnhj

Hii baptiste,

I tried you example, but I get still one plot and not as you said two plots
...
I use R for Windows, could this have something to do with this problem ?



baptiste auguie-2 wrote:
> 
> 
> On 23 Mar 2009, at 11:52, johnhj wrote:
> 
>>
>> I have still the same problem... As you said I tried with  
>> par(mfrow=c(2,1))
>> and par(mfrow=c(1,2)) but without success. Could R compiler be the  
>> problem ?
>>
> 
> Why not? But may I suggest you try first the following example I sent  
> you yesterday,
> 
> 
>> png()
>>
>> par(mfrow=c(2,1) )
>>
>> plot(1:10)
>>
>> plot(1:10)
>>
>> dev.off()
> 
> 
> Does this not produce two plots on a single page?
> 
> 
> baptiste
> 
> 
> PS: please do read the posting guide
>>
>> R-help@r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide
>> http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
> 
> 
> 
>>
>> Wills, Kellie wrote:
>>>
>>> par(mfrow=c(1,1)) will give you just one panel.  Try  
>>> par(mfrow=c(2,1)) or
>>> par(mfrow=c(1,2)).
>>>
>>>
>>> -Original Message-
>>> From: r-help-boun...@r-project.org on behalf of johnhj
>>> Sent: Sun 3/22/2009 10:50 AM
>>> To: r-help@r-project.org
>>> Subject: [R]  Problems with combining plots
>>>
>>>
>>> Hii,
>>>
>>> I will combine some plots. Like this example here
>>> http://www.statmethods.net/advgraphs/layout.html I tired to do it  
>>> for 2
>>> plots but without success.
>>>
>>> Here is my code:
>>>
>>> test<-read.table(file="D:/file.txt")
>>> space<-read.table(file="D:/space.txt")
>>>
>>> space$gruppe <- 502*rep(1:6, each=7)
>>> x<- c(test$V1)
>>> y<- c(test$V2)
>>>
>>> par(mfrow=c(1,1))
>>>
>>> png(filename = "D:/example.png", width = 640, height = 480,  
>>> pointsize =
>>> 12,
>>> bg = "white",  res = NA)
>>>
>>> boxplot(V2 ~ gruppe , data = space , col = "lightgray",boxwex=0.2)
>>> plot(panel.first=grid(ny=NULL,nx=NULL),x,y, xlab = "Zeit(sec)", ylab
>>> ="Datenrate(MBit(sec))",ylim=c(0,40), col ="purple", type ="l", main
>>> ="combined plots",lwd=2)
>>>
>>> dev.off()
>>>
>>>
>>> What is the mistake in my code ?
>>>
>>>
>>> -- 
>>> View this message in context:
>>> http://www.nabble.com/Problems-with-combining-plots-tp22646692p22646692.html
>>> Sent from the R help mailing list archive at Nabble.com.
>>>
>>> __
>>> R-help@r-project.org mailing list
>>> https://stat.ethz.ch/mailman/listinfo/r-help
>>> PLEASE do read the posting guide
>>> http://www.R-project.org/posting-guide.html
>>> and provide commented, minimal, self-contained, reproducible code.
>>>
>>> __
>>> R-help@r-project.org mailing list
>>> https://stat.ethz.ch/mailman/listinfo/r-help
>>> PLEASE do read the posting guide
>>> http://www.R-project.org/posting-guide.html
>>> and provide commented, minimal, self-contained, reproducible code.
>>>
>>>
>>
>> -- 
>> View this message in context:
>> http://www.nabble.com/Problems-with-combining-plots-tp22646692p22658673.html
>> Sent from the R help mailing list archive at Nabble.com.
>>
>> __
>> R-help@r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide
>> http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
> 
> _
> 
> Baptiste Auguié
> 
> School of Physics
> University of Exeter
> Stocker Road,
> Exeter, Devon,
> EX4 4QL, UK
> 
> Phone: +44 1392 264187
> 
> http://newton.ex.ac.uk/research/emag
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 
> 

-- 
View this message in context: 
http://www.nabble.com/Problems-with-combining-plots-tp22646692p22662039.html
Sent from the R help mailing list archive at Nabble.com.

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


Re: [R] Copying a subset of a matrix

2009-03-23 Thread akintayo holder
Thanks. That worked.

On Mon, Mar 23, 2009 at 2:17 PM, Phil Spector wrote:

> Suppose your matrix is called "x".  Then
>
>x[,seq(1,ncol(x),n)]
>
> will give you what you want.
>
>   - Phil Spector
> Statistical Computing Facility
> Department of Statistics
> UC Berkeley
> spec...@stat.berkeley.edu
>
>
>
> On Mon, 23 Mar 2009, akintayo holder wrote:
>
>  Hi,
>> I have a matrix and I want to create a matrix that includes every nth
>> column
>> of the original matrix. Does anyone know how I can go about doing that ?
>> Or
>> if you have an idea how I can do the same thing with a data frame ? If you
>> can just point me towards the appropriate function I would appreciate it.
>> Thanks
>>
>>[[alternative HTML version deleted]]
>>
>> __
>> R-help@r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide
>> http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>>
>>

[[alternative HTML version deleted]]

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


Re: [R] Basic regression output question

2009-03-23 Thread Kingsford Jones
On Mon, Mar 23, 2009 at 11:04 AM, GRANT Lewis  wrote:

[snip]

> factors<-matrix(runif(400),nrow=40)
> returns1<-matrix(runif(40),nrow=1)
> returns2<-matrix(runif(2000),nrow=50)
> coef(summary(lm(t(returns1)~factors)))[1,4]

[snip]

> (coef(summary(lm(t(returns2)~factors)))[50])[1,4]
>
> Error in `[.default`((coef(summary(lm(t(returns2) ~ factors)))[50]), 1,
> :
>  incorrect number of dimensions


Use of the 'str' function will reveal you are indexing a list element,
so try [[50]] rather than [50]

hth,

Kingsford Jones



>
>
>
> **
> Hermes Fund Managers Limited
> Registered in England No. 1661776, Lloyds Chambers, 1 Portsoken Street, 
> London E1 8HZ
>
> *** Please read the Hermes email disclaimer at 
> http://www.hermes.co.uk/email_terms.htm before acting on this email or 
> opening any attachment ***
>
> The contents of this email are confidential.  If you have received this 
> message in error, please delete it immediately and contact the sender 
> directly or the Hermes IT Helpdesk on +44(0)20 7680 2117.  Any reliance on, 
> use, disclosure, dissemination, distribution or copying of this email is 
> unauthorised and strictly prohibited.
>
> This message has been checked for viruses but the recipient is strongly 
> advised to rescan the message before opening any attachments or attached 
> executable files.  Hermes do not accept any liability for any damage 
> sustained as a result of a virus introduced by this email or any attachment.
>
>
> **
>
>
> __
> This email has been scanned by the MessageLabs Email Security System.
> For more information please visit http://www.messagelabs.com/email
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

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


Re: [R] How to set up a function for "Central Limit Theorem"

2009-03-23 Thread Charles Annis, P.E.
Ivan:

While you're figuring out how to execute the CLT in R you may find my
automated examples informative.

http://StatisticalEngineering.com/central_limit_theorem.htm

Charles Annis, P.E.

charles.an...@statisticalengineering.com
phone: 561-352-9699
eFax:  614-455-3265
http://www.StatisticalEngineering.com
 

-Original Message-
From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On
Behalf Of pfc_ivan
Sent: Monday, March 23, 2009 1:31 PM
To: r-help@r-project.org
Subject: [R] How to set up a function for "Central Limit Theorem"


Hello guys, I am stuck here:

How do I make 1000 samples of n = 10 observations from an Exponential
distribution and then compute the mean for all those 1000 samples? 

Basically I need to prove the Central Limit theorem, which states:

http://www.nabble.com/file/p22664113/d175f06cbf200bd52a2c27a2e56dc594.png 

Where the Sn is sum of random variables, n we have from the question, mu is
mean and (sigma)^2 is variance. 

I am having trouble setting up the function to do this. 

Any help apreciated!

-- 
View this message in context:
http://www.nabble.com/How-to-set-up-a-function-for-%22Central-Limit-Theorem%
22-tp22664113p22664113.html
Sent from the R help mailing list archive at Nabble.com.

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

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


Re: [R] Extracting SD of random effects from lme object

2009-03-23 Thread Kingsford Jones
On Mon, Mar 23, 2009 at 11:26 AM, Ben Domingue  wrote:
> Hello,
> How do I get the standard deviations for the random effects out of the
> lme object?  I feel like there's probably a simple way of doing this,
> but I can't see it.  Using the first example from the documentation:
>
>> fm1 <- lme(distance ~ age, data = Orthodont) # random is ~ age
>> fm1
> Linear mixed-effects model fit by REML
>  Data: Orthodont
>  Log-restricted-likelihood: -221.3183
>  Fixed: distance ~ age
> (Intercept)         age
>  16.761   0.6601852
>
> Random effects:
>  Formula: ~age | Subject
>  Structure: General positive-definite
>            StdDev    Corr
> (Intercept) 2.3270339 (Intr)
> age         0.2264276 -0.609
> Residual    1.3100399
>
> Number of Observations: 108
> Number of Groups: 27
>
> I want to extract the column vector (2.3270339, 0.2264276,
> 1.3100399)'.  Any thoughts?

To get the covariance matrix of the random effects:

> getVarCov(fm1)
Random effects variance covariance matrix
(Intercept)   age
(Intercept) 5.41510 -0.321060
age-0.32106  0.051269
  Standard Deviations: 2.327 0.22643


One way to extract the standard deviations shown by the print method above is:

> diag(sqrt(getVarCov(fm1)))
(Intercept) age
  2.3270339   0.2264276
Warning message:
In sqrt(getVarCov(fm1)) : NaNs produced


And to get the estimate of the error standard deviation:

> fm1$sigma
[1] 1.31004


hth,

Kingsford Jones

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

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


Re: [R] Extracting SD of random effects from lme object

2009-03-23 Thread Ben Domingue
On Mon, Mar 23, 2009 at 1:18 PM, Kingsford Jones
 wrote:
> On Mon, Mar 23, 2009 at 11:26 AM, Ben Domingue  wrote:
>> Hello,
>> How do I get the standard deviations for the random effects out of the
>> lme object?  I feel like there's probably a simple way of doing this,
>> but I can't see it.  Using the first example from the documentation:
>>
>>> fm1 <- lme(distance ~ age, data = Orthodont) # random is ~ age
>>> fm1
>> Linear mixed-effects model fit by REML
>>  Data: Orthodont
>>  Log-restricted-likelihood: -221.3183
>>  Fixed: distance ~ age
>> (Intercept)         age
>>  16.761   0.6601852
>>
>> Random effects:
>>  Formula: ~age | Subject
>>  Structure: General positive-definite
>>            StdDev    Corr
>> (Intercept) 2.3270339 (Intr)
>> age         0.2264276 -0.609
>> Residual    1.3100399
>>
>> Number of Observations: 108
>> Number of Groups: 27
>>
>> I want to extract the column vector (2.3270339, 0.2264276,
>> 1.3100399)'.  Any thoughts?
>
> To get the covariance matrix of the random effects:
>
>> getVarCov(fm1)
> Random effects variance covariance matrix
>            (Intercept)       age
> (Intercept)     5.41510 -0.321060
> age            -0.32106  0.051269
>  Standard Deviations: 2.327 0.22643
>
>
> One way to extract the standard deviations shown by the print method above is:
>
>> diag(sqrt(getVarCov(fm1)))
> (Intercept)         age
>  2.3270339   0.2264276
> Warning message:
> In sqrt(getVarCov(fm1)) : NaNs produced
>
>
> And to get the estimate of the error standard deviation:
>
>> fm1$sigma
> [1] 1.31004
>
>
> hth,
>
That's what I needed.  Thanks.


Ben


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

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


[R] specifying repeated measures model in lmer

2009-03-23 Thread Lawrence Hanser
Dear Colleagues,
I have what Roger Kirk (Experimental Design: Procedures for the Behavioral
Sciences, 1968) refers to as a randomized block factorial design.  The anova
table would look like this:

  df
A 3
Subj/A  103 (error term for A)
B   23
A*B69
B*Subj/A 2369 (error term for B and A*B)

Subjects are nested within A and give a response for each B.  If y is the
dependent variable, is this the correct lmer specification for the above,
where ID is the variable name for Subj:

lmer(y ~ A + B + A*B + (A|ID))

Am I barking up the right tree?  I can also fit:

aov(y ~ A + B + A*B  + ID)
then I have to do some hand calculations to use ID as the error term for A.
 The residual (really B*ID) is the correct error term for B and A*B.

Thanks,

Larry

[[alternative HTML version deleted]]

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


Re: [R] How to set up a function for "Central Limit Theorem"

2009-03-23 Thread Greg Snow
There is also the clt.examp function in the TeachingDemos package.

-- 
Gregory (Greg) L. Snow Ph.D.
Statistical Data Center
Intermountain Healthcare
greg.s...@imail.org
801.408.8111


> -Original Message-
> From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-
> project.org] On Behalf Of Charles Annis, P.E.
> Sent: Monday, March 23, 2009 1:14 PM
> To: 'pfc_ivan'; r-help@r-project.org
> Subject: Re: [R] How to set up a function for "Central Limit Theorem"
> 
> Ivan:
> 
> While you're figuring out how to execute the CLT in R you may find my
> automated examples informative.
> 
> http://StatisticalEngineering.com/central_limit_theorem.htm
> 
> Charles Annis, P.E.
> 
> charles.an...@statisticalengineering.com
> phone: 561-352-9699
> eFax:  614-455-3265
> http://www.StatisticalEngineering.com
> 
> 
> -Original Message-
> From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-
> project.org] On
> Behalf Of pfc_ivan
> Sent: Monday, March 23, 2009 1:31 PM
> To: r-help@r-project.org
> Subject: [R] How to set up a function for "Central Limit Theorem"
> 
> 
> Hello guys, I am stuck here:
> 
> How do I make 1000 samples of n = 10 observations from an Exponential
> distribution and then compute the mean for all those 1000 samples?
> 
> Basically I need to prove the Central Limit theorem, which states:
> 
> http://www.nabble.com/file/p22664113/d175f06cbf200bd52a2c27a2e56dc594.p
> ng
> 
> Where the Sn is sum of random variables, n we have from the question,
> mu is
> mean and (sigma)^2 is variance.
> 
> I am having trouble setting up the function to do this.
> 
> Any help apreciated!
> 
> --
> View this message in context:
> http://www.nabble.com/How-to-set-up-a-function-for-%22Central-Limit-
> Theorem%
> 22-tp22664113p22664113.html
> Sent from the R help mailing list archive at Nabble.com.
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-
> guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-
> guide.html
> and provide commented, minimal, self-contained, reproducible code.

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


Re: [R] Extracting SD of random effects from lme object

2009-03-23 Thread Douglas Bates
On Mon, Mar 23, 2009 at 2:18 PM, Kingsford Jones
 wrote:
> On Mon, Mar 23, 2009 at 11:26 AM, Ben Domingue  wrote:
>> Hello,
>> How do I get the standard deviations for the random effects out of the
>> lme object?  I feel like there's probably a simple way of doing this,
>> but I can't see it.  Using the first example from the documentation:
>>
>>> fm1 <- lme(distance ~ age, data = Orthodont) # random is ~ age
>>> fm1
>> Linear mixed-effects model fit by REML
>>  Data: Orthodont
>>  Log-restricted-likelihood: -221.3183
>>  Fixed: distance ~ age
>> (Intercept)         age
>>  16.761   0.6601852
>>
>> Random effects:
>>  Formula: ~age | Subject
>>  Structure: General positive-definite
>>            StdDev    Corr
>> (Intercept) 2.3270339 (Intr)
>> age         0.2264276 -0.609
>> Residual    1.3100399
>>
>> Number of Observations: 108
>> Number of Groups: 27
>>
>> I want to extract the column vector (2.3270339, 0.2264276,
>> 1.3100399)'.  Any thoughts?
>
> To get the covariance matrix of the random effects:
>
>> getVarCov(fm1)
> Random effects variance covariance matrix
>            (Intercept)       age
> (Intercept)     5.41510 -0.321060
> age            -0.32106  0.051269
>  Standard Deviations: 2.327 0.22643
>
>
> One way to extract the standard deviations shown by the print method above is:
>
>> diag(sqrt(getVarCov(fm1)))
> (Intercept)         age
>  2.3270339   0.2264276
> Warning message:
> In sqrt(getVarCov(fm1)) : NaNs produced

You can avoid the warning message if you extract the diagonal first
then take the square root.

sqrt(diag(getVarCov(fm1)))

>
>
> And to get the estimate of the error standard deviation:
>
>> fm1$sigma
> [1] 1.31004
>
>
> hth,
>
> Kingsford Jones
>
>> Thanks,
>> Ben
>>
>> __
>> R-help@r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>>
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

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


Re: [R] specifying repeated measures model in lmer

2009-03-23 Thread Douglas Bates
On Mon, Mar 23, 2009 at 2:35 PM, Lawrence Hanser  wrote:
> Dear Colleagues,
> I have what Roger Kirk (Experimental Design: Procedures for the Behavioral
> Sciences, 1968) refers to as a randomized block factorial design.  The anova
> table would look like this:
>
>                      df
> A                     3
> Subj/A          103 (error term for A)
> B                   23
> A*B                69
> B*Subj/A     2369 (error term for B and A*B)

> Subjects are nested within A and give a response for each B.  If y is the
> dependent variable, is this the correct lmer specification for the above,
> where ID is the variable name for Subj:

> lmer(y ~ A + B + A*B + (A|ID))

If, as you say, subjects are nested within levels of A, then I don't
think you want a random effects term of the form (A | ID).  I
understand what you say to mean that each subject is exposed to one
and only one level of factor A so trying to fit a random effect for
the levels of A within each subject doesn't make sense.

Trying to understand model specifications for lmer according to the
degrees of freedom for each term is probably not the best approach.

> Am I barking up the right tree?  I can also fit:
>
> aov(y ~ A + B + A*B  + ID)
> then I have to do some hand calculations to use ID as the error term for A.
>  The residual (really B*ID) is the correct error term for B and A*B.
>
> Thanks,
>
> Larry
>
>        [[alternative HTML version deleted]]
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

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


Re: [R] Difference between gam() and loess().

2009-03-23 Thread Rolf Turner


On 21/03/2009, at 3:19 AM, Ravi Varadhan wrote:




I also tried a number of other things including changing the  
"family", and parameters in
"loess.control", but to no avail.  I looked at the Fortran codes  
from both loess and gam.
They are daunting, to say the least. They are dense, and there are  
absolutely no comments
whatsoever.  But one thing is clear - they are using different  
Fortran codes.


	Thanks for doing this digging.  It would on this basis be not  
unreasonable

to expect there to be numerical differences in the result.


So, the best bet might be to get Trevor Hastie or Bill Cleveland to  
help you out.


But, before that:  why is this an issue, Rolf?  Is it important  
that these two results

be identical?


	It just makes me nervous when two procedures which I believe to be  
doing the
	same thing give answers which are not identical.  Such a phenomenon  
makes me
	wonder if there is something which I am not understanding, and if  
thereby my

lack of understanding might lead to my making serious errors.

	A lack of understanding is an all-to-frequent occurrence in my  
circumstances
	(i.e. being thick as two short planks) and I therefore need to be  
constantly

on my guard against such a lack.

	In the current setting I think I will desist from pursuing the issue  
any further
	and just assume that different coding accounts for the different  
results.  And

take a valium.

cheers,

Rolf Turner

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

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


[R] aggregate() problem on last patched version

2009-03-23 Thread Kenneth Roy Cabrera Torres
Hi R users and developers:

Several hours ago I post a problem with the aggregate() function
using the last patch version of R.
(R version 2.8.1 Patched (2009-03-18 r48193))

I do compile again the 2.8.1 plain version (not patched) and
now it works again, both the examples and  my old scripts that
use this function.

What are the changes on the patched version that produces this
behavior on this particular function?

Thank you for your attention.

Kenneth

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


[R] Graphic with several curves

2009-03-23 Thread CE.KA

Hi R users,

Imagine the folowing data frame
  19901991   1992
1   5  20  6
2   15  1  11
3   3   14 22
4   20  8  55
5   10  3  14

Is there a way to build a graphic in which:
- 1 curve represents the observation 1 during 1990 & 1992
- 1 curve represents the observation 2 during 1990 & 1992
- 1 curve represents the observation 3 during 1990 & 1992
- 1 curve represents the observation 4 during 1990 & 1992
- 1 curve represents the observation 5 during 1990 & 1992

Sincerely yours
-- 
View this message in context: 
http://www.nabble.com/Graphic-with-several-curves-tp22669675p22669675.html
Sent from the R help mailing list archive at Nabble.com.

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


[R] read in large data file (tsv) with inline filter?

2009-03-23 Thread David Reiss
I have a very large tab-delimited file, too big to store in memory via
readLines() or read.delim(). Turns out I only need a few hundred of those
lines to be read in. If it were not so large, I could read the entire file
in and "grep" the lines I need. For such a large file; many calls to
read.delim() with incrementing "skip" and "nrows" parameters, followed by
grep() calls is very slow. I am aware of possibilities via SQLite; I would
prefer to not use that in this case.

My question is...Is there a function for efficiently reading in a file along
the lines of read.delim(), which allows me to specify a filter (via grep or
something else) that tells the function to only read in certain lines that
match?

If not, I would *love* to see a "filter" parameter added as an option to
read.delim() and/or readLines().

thanks for any pointers.

--David

[[alternative HTML version deleted]]

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


[R] Merging rows in dataframes

2009-03-23 Thread Schraga Schwartz

Hello,

I have a dataframe with 40 columns and around 450,000 rows. The first  
column in each row is a factor id and the remaining are numeric. Some  
rows have the same ids. What I want to do is to merge each set of rows  
sharing the same ids (id set) into one single row (summarizing row)  
with that id. To create the summarizing row, I'd like to apply a  
different function on each of the original columns in the id set. Some  
columns within the summarizing row will equal the mean of the columns  
in the id set, others will equal the minimum, others the maximum.


To do this, I tried using the by() function. However, this was  
extremely slow (it ran for more than two hours before I stopped it).  
Also, it used up all of 16 GB of memory on my machine. Is there any  
more efficient function, both in terms of time and memory, to do this  
sort of thing?


Thank you very much,
Schraga Schwartz

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


[R] [Fwd: Re: Graphic with several curves]

2009-03-23 Thread Kenneth Roy Cabrera Torres
- Mensaje reenviado 
> De: Kenneth Roy Cabrera Torres 
> Para: CE.KA 
> Asunto: Re: [R] Graphic with several curves
> Fecha: Mon, 23 Mar 2009 17:02:48 -0500
> 
> try matplot(),
> 
> but you have to transpose the matrix
> 
> with option type="l".
> 
> Example:
> 
> data1<-data.frame(X1990=c(5,15,3,20,10),X1991=c(20,1,14,8,3),X1992=c(6,11,22,55,14))
> data2<-t(data1)
> matplot(as.integer(substring(row.names(data2),2)),data2,type="l")
> 
> Data frames names cannot start with a number. R change it to an X,
> then when you transpose I obtain the substring with the number that I
> need.
> 
> Hope, this helps (HTH)
> 
> Kenneth
> 
> 
> 
> 
> 
> El lun, 23-03-2009 a las 14:33 -0700, CE.KA escribió:
> > Hi R users,
> > 
> > Imagine the folowing data frame
> >   19901991   1992
> > 1   5  20  6
> > 2   15  1  11
> > 3   3   14 22
> > 4   20  8  55
> > 5   10  3  14
> > 
> > Is there a way to build a graphic in which:
> > - 1 curve represents the observation 1 during 1990 & 1992
> > - 1 curve represents the observation 2 during 1990 & 1992
> > - 1 curve represents the observation 3 during 1990 & 1992
> > - 1 curve represents the observation 4 during 1990 & 1992
> > - 1 curve represents the observation 5 during 1990 & 1992
> > 
> > Sincerely yours

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


[R] Confusion regarding environments invoked by "source" command

2009-03-23 Thread Dennis Fisher

Colleagues,

R version 2.8.1 in OS X

Within a function (which is already within a function), I am sourcing  
a file.  The syntax of the command is something like (this is just an  
example; the actual code is much more complicated):

BIGFUNCTION <- function()
{
DATAFRAME   <- [some commands to create a dataframe]
MYFUNCTION(DATAFRAME)
}
MYFUNCTION  <- function(DATAFRAME)
{
print(ls())
exists("DATAFRAME")
source("myfile", local=T)
}

The file "myfile" contains the following:
print(DATAFRAME)

When I execute BIGFUNCTION from the command line, R reports:
Error in print(DATAFRAME) : object "DATAFRAME" not found
	Calls: BIGFUNCTION ... source -> eval.with.vis -> eval.with.vis ->  
print -> head

Execution halted
However, the "print(ls())" command has just listed the existence of  
DATAFRAME and "exists(DATAFRAME) return TRUE
So, I am at a loss as to why ls() and "exists" see the object but  
print(DATAFRAME) command does not.


I tried "local=F" in the source command to no avail.
I also replaced "print(DATAFRAME)" with "print(str(DATAFRAME))" with  
the same outcome.
I assume that this is an error as to how I am managing environments  
but I cannot figure out a solution.  Any help would be appreciated.


Dennis



Dennis Fisher MD
P < (The "P Less Than" Company)
Phone: 1-866-PLessThan (1-866-753-7784)
Fax: 1-415-564-2220
www.PLessThan.com

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


Re: [R] read in large data file (tsv) with inline filter?

2009-03-23 Thread Dylan Beaudette
On Monday 23 March 2009, David Reiss wrote:
> I have a very large tab-delimited file, too big to store in memory via
> readLines() or read.delim(). Turns out I only need a few hundred of those
> lines to be read in. If it were not so large, I could read the entire file
> in and "grep" the lines I need. For such a large file; many calls to
> read.delim() with incrementing "skip" and "nrows" parameters, followed by
> grep() calls is very slow. I am aware of possibilities via SQLite; I would
> prefer to not use that in this case.
>
> My question is...Is there a function for efficiently reading in a file
> along the lines of read.delim(), which allows me to specify a filter (via
> grep or something else) that tells the function to only read in certain
> lines that match?
>
> If not, I would *love* to see a "filter" parameter added as an option to
> read.delim() and/or readLines().
>
> thanks for any pointers.
>
> --David


How about pre-filtering before loading the data into R:

grep -E 'your pattern here' your_file_here > your_filtered_file

alternatively if you need to search in fields, see 'awk', and 'cut', or if you 
need to delete things see 'tr'.

These tools come with any unix-like OS, and you can probably get them on 
windows without much effort.


Cheers,
Dylan


-- 
Dylan Beaudette
Soil Resource Laboratory
http://casoilresource.lawr.ucdavis.edu/
University of California at Davis
530.754.7341

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


Re: [R] Graphic with several curves

2009-03-23 Thread Jorge Ivan Velez
Dear CE.KA,
If 'x' is your data, you could do something like the following:

matplot(as.factor(c('1990','1992')),t(x[,-2]),type='l',
   xaxt='n',xlab='Year',ylab="Some label here",lty=1)
axis(1,c(1990,"","","",1992),c(1990,"","","",1992))
legend('topleft',paste('Observation',1:5),col=1:5,lty=1,cex=0.8)

See ?matplot, ?axis and ?legend for more details.

HTH,

Jorge


On Mon, Mar 23, 2009 at 5:33 PM, CE.KA  wrote:

>
> Hi R users,
>
> Imagine the folowing data frame
>  19901991   1992
> 1   5  20  6
> 2   15  1  11
> 3   3   14 22
> 4   20  8  55
> 5   10  3  14
>
> Is there a way to build a graphic in which:
> - 1 curve represents the observation 1 during 1990 & 1992
> - 1 curve represents the observation 2 during 1990 & 1992
> - 1 curve represents the observation 3 during 1990 & 1992
> - 1 curve represents the observation 4 during 1990 & 1992
> - 1 curve represents the observation 5 during 1990 & 1992
>
> Sincerely yours
> --
> View this message in context:
> http://www.nabble.com/Graphic-with-several-curves-tp22669675p22669675.html
> Sent from the R help mailing list archive at Nabble.com.
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

[[alternative HTML version deleted]]

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


Re: [R] Accuracy of R and other platforms

2009-03-23 Thread John Maindonald
These comparisons are very simplistic.  In most contexts, it would  
make much better sense to measure "accuracy" in standard error units,  
rather than in number of digits.

There doubtless are specialist applications where the 15th digit (or  
even the 10th!) are important.  But the check of accuracy really has  
then to be tuned to that application.  Results in cases where the  
"accuracy" beyond (optimistically) the sixth digit is of no  
consequence are unlikely to give much clue on performance in such  
specialist applications.

John Maindonald email: john.maindon...@anu.edu.au
phone : +61 2 (6125)3473fax  : +61 2(6125)5549
Centre for Mathematics & Its Applications, Room 1194,
John Dedman Mathematical Sciences Building (Building 27)
Australian National University, Canberra ACT 0200.


On 23/03/2009, at 10:00 PM, r-help-requ...@r-project.org wrote:

> From: Karl Ove Hufthammer 
> Date: 23 March 2009 2:58:23 AM
> To: r-h...@stat.math.ethz.ch
> Subject: Re: [R] Accuracy of R and other platforms
>
>
> Alejandro C. Frery:
>
>> @ARTICLE{AlmironSilvaMM:2009,
>> author = {Almiron, M. and Almeida, E. S. and Miranda, M.},
>> title = {The Reliability of Statistical Functions in Four Software
>> Packages Freely used in Numerical Computation},
>> journal = {Brazilian Journal of Probability and Statistics},
>> year = {in press},
>> volume = {Special Issue on Statistical Image and Signal Processing},
>> url = {http://www.imstat.org/bjps/}}
>>
>> is freely available under the "Future Papers" link. It makes a nice
>> comparison of the numerical properties of R, Ox, Octave and Python.
>
> Thanks for posting this. I’m happy to see that the results for R were
> generally excellent, and almost always better than for the three other
> software packages.
>
> But there were a few cases where R did not turn out to be the winner.
> Rather surprising that Ox was better than R for computing the
> autocorrelation coefficient for two of the datasets, given its  
> terrible
> results for the standard deviation. Anybody have any ideas why?
>
> -- 
> Karl Ove Hufthammer


[[alternative HTML version deleted]]

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


Re: [R] Confusion regarding environments invoked by "source" command

2009-03-23 Thread Duncan Murdoch

On 23/03/2009 6:06 PM, Dennis Fisher wrote:

Colleagues,

R version 2.8.1 in OS X

Within a function (which is already within a function), I am sourcing  
a file.  The syntax of the command is something like (this is just an  
example; the actual code is much more complicated):


This code works just as expected, once you make it syntactically 
correct.  I suspect you have something else wrong in the real case. 
Here's my log:


> BIGFUNCTION<- function()
+ {
+ DATAFRAME<- data.frame(a=1:3,b=4:6)
+ MYFUNCTION(DATAFRAME)
+ }
>
> MYFUNCTION<- function(DATAFRAME)
+ {
+ print(ls())
+ exists("DATAFRAME")
+ source("myfile.R", local=T)
+ }
> BIGFUNCTION()
[1] "DATAFRAME"
  a b
1 1 4
2 2 5
3 3 6

Duncan Murdoch



BIGFUNCTION <- function()
{
DATAFRAME   <- [some commands to create a dataframe]
MYFUNCTION(DATAFRAME)
}
MYFUNCTION  <- function(DATAFRAME)
{
print(ls())
exists("DATAFRAME")
source("myfile", local=T)
}

The file "myfile" contains the following:
print(DATAFRAME)

When I execute BIGFUNCTION from the command line, R reports:
Error in print(DATAFRAME) : object "DATAFRAME" not found
	Calls: BIGFUNCTION ... source -> eval.with.vis -> eval.with.vis ->  
print -> head

Execution halted
However, the "print(ls())" command has just listed the existence of  
DATAFRAME and "exists(DATAFRAME) return TRUE
So, I am at a loss as to why ls() and "exists" see the object but  
print(DATAFRAME) command does not.


I tried "local=F" in the source command to no avail.
I also replaced "print(DATAFRAME)" with "print(str(DATAFRAME))" with  
the same outcome.
I assume that this is an error as to how I am managing environments  
but I cannot figure out a solution.  Any help would be appreciated.


Dennis



Dennis Fisher MD
P < (The "P Less Than" Company)
Phone: 1-866-PLessThan (1-866-753-7784)
Fax: 1-415-564-2220
www.PLessThan.com

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


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


[R] Plot Means Line with Standard Deviation as "Whiskers"

2009-03-23 Thread Rodrigo Aluizio
Hi list members.

I’ll try to plot the abundance means of nine transects as lines, with five
points on each transect (A to I). I will also need to plot for each point,
it’s standard deviation (once each point will have tree replicates) as
whiskers. Another problem will be that all the points should exist for each
transect, but we know that six of the nine transects will have blanks (A1,
A2, A3, -, A5). So the line have to be interrupted and can’t induce an
inexistent value by link A3 to A5 (an example).

 

Well, unfortunately I still don’t have the data, I’m just anticipating a
future situation, and I can’t find any clear suggestion by “googling” the
issue, and I guess the help files may have something but I was unable to
find by myself.

 

Thanks in advance. Any help will be appreciated.

 

-

MSc.   Rodrigo Aluizio

Centro de Estudos do Mar/UFPR
Laboratório de Micropaleontologia
Avenida Beira Mar s/n - CEP 83255-000
Pontal do Paraná - PR - Brasil
Fone: (0**41) 3455-3620 ramal 217
Fax: (0**41) 3455-3623




[[alternative HTML version deleted]]

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


Re: [R] .Internal

2009-03-23 Thread rkevinburton
Sorry to be so dense but the article that you suggest does not give any 
information on how the arguments are packed up. I look at the call:

val <- .Internal(fmin(function(arg) -f(arg, ...), lower, upper, tol))

and then with the help of this article I find do_fmin in optimize.c:

SEXP attribute_hidden do_fmin(SEXP call, SEXP op, SEXP args, SEXP rho)

Again there doesn't seem to be any coorespondance between lower, upper, tol and 
the arguments to do_fmin. So I am missing a step.

Thank you.

Kevin


 Berwin A Turlach  wrote: 
> G'day Kevin,
> 
> On Wed, 18 Mar 2009 21:46:51 -0700
>  wrote:
> 
> > I was trying to find source for optimize and I ran across
> > 
> > function (f, interval, ..., lower = min(interval), upper =
> > max(interval), maximum = FALSE, tol = .Machine$double.eps^0.25) 
> > {
> > if (maximum) {
> > val <- .Internal(fmin(function(arg) -f(arg, ...), lower, 
> > upper, tol))
> > list(maximum = val, objective = f(val, ...))
> > }
> > else {
> > val <- .Internal(fmin(function(arg) f(arg, ...), lower, 
> > upper, tol))
> > list(minimum = val, objective = f(val, ...))
> > }
> > }
> > 
> > Then I did a search for fmin and i came up with:
> > 
> > /* fmin(f, xmin, xmax tol) */
> > SEXP attribute_hidden do_fmin(SEXP call, SEXP op, SEXP args, SEXP rho)
> > 
> > 
> > So my question is where do I find the intermediary step between 
> > 
> > .Internal(fmin(function(arg) f(arg, ...), lower,  upper, tol))
> > 
> > and
> > 
> > SEXP attribute_hidden do_fmin(SEXP call, SEXP op, SEXP args, SEXP rho)
> 
> @Article{Rnews:Ligges:2006,
>   author   = {Uwe Ligges},
>   title= {{R} {H}elp {D}esk: {Accessing} the Sources},
>   journal  = {R News},
>   year = 2006,
>   volume   = 6,
>   number   = 4,
>   pages= {43--45},
>   month= {October},
>   url  = http,
>   pdf= Rnews2006-4
> }
> 
> http://CRAN.R-project.org/doc/Rnews/Rnews_2006-4.pdf
> 
> > The number of arguments doesn't match up. I am guessing that lower
> > and upper somehow get merged into the args. And rho is 'tol'. Right?
> 
> Unlikely.  In "Writing R Extensions" (and the functions I looked up),
> 'rho' usually denotes an environment that is used to evaluate
> expressions in.  Typically (i.e. in cases that I had need to look at),
> all arguments are rolled into the SEXP arg for internal functions.
> 
> HTH.
> 
> Cheers,
> 
>   Berwin

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


Re: [R] Merging rows in dataframes

2009-03-23 Thread Gabor Grothendieck
Using sqldf you only need two statements, infile <- file(...) and
DF <- sqldf("select min(a), max(b), mean(c), ... from infile group by id").
The file statement identifies the filename and the second reads it
into sqlite (without
going through R), summarizes it and then reads the summarized version
into R.  You may also need to provide info on its format if its not in the
default format.  See example 4a on home page and the other examples
there:
http://sqldf.googlecode.com


On Mon, Mar 23, 2009 at 5:58 PM, Schraga Schwartz
 wrote:
> Hello,
>
> I have a dataframe with 40 columns and around 450,000 rows. The first column
> in each row is a factor id and the remaining are numeric. Some rows have the
> same ids. What I want to do is to merge each set of rows sharing the same
> ids (id set) into one single row (summarizing row) with that id. To create
> the summarizing row, I'd like to apply a different function on each of the
> original columns in the id set. Some columns within the summarizing row will
> equal the mean of the columns in the id set, others will equal the minimum,
> others the maximum.
>
> To do this, I tried using the by() function. However, this was extremely
> slow (it ran for more than two hours before I stopped it). Also, it used up
> all of 16 GB of memory on my machine. Is there any more efficient function,
> both in terms of time and memory, to do this sort of thing?
>
> Thank you very much,
> Schraga Schwartz
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

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


[R] Plotting with plot()

2009-03-23 Thread Satindra Chakravorty
I am trying to use the plot() function in R but can't seem to generate any
plots.  I am running R version 2.7.1 on a server machine to which I am
remotely connected on a Linux platform.

I typed
> plot(trees$volume)
as a test and got back the ">" prompt with no plot in my current window or
in any other window.

I read in another post on this list about opening a graphics device with
x11().  I tried that and got an ERROR message:
Error in X11
   unable to start device X11
In x11() : unable to open connection to X11 display ' '.

Please let me know what I should do.

Thanks for your help.
Satindra.

[[alternative HTML version deleted]]

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


Re: [R] Plot Means Line with Standard Deviation as "Whiskers"

2009-03-23 Thread Peter Alspach
Tena koe Rodrigo

You could always make up some data and then show us what you have tried to do.  
I would guess you need to check out:

plot # to do the basic plot
lines # to add lines to the plot
points # to add points to the plot
arrows # can be used to give the whiskers
apply # to get the means and standard deviations

HTH 

Peter Alspach

> -Original Message-
> From: r-help-boun...@r-project.org 
> [mailto:r-help-boun...@r-project.org] On Behalf Of Rodrigo Aluizio
> Sent: Tuesday, 24 March 2009 11:44 a.m.
> To: R Help
> Subject: [R] Plot Means Line with Standard Deviation as "Whiskers"
> 
> Hi list members.
> 
> I'll try to plot the abundance means of nine transects as 
> lines, with five points on each transect (A to I). I will 
> also need to plot for each point, it's standard deviation 
> (once each point will have tree replicates) as whiskers. 
> Another problem will be that all the points should exist for 
> each transect, but we know that six of the nine transects 
> will have blanks (A1, A2, A3, -, A5). So the line have to be 
> interrupted and can't induce an inexistent value by link A3 
> to A5 (an example).
> 
>  
> 
> Well, unfortunately I still don't have the data, I'm just 
> anticipating a future situation, and I can't find any clear 
> suggestion by "googling" the issue, and I guess the help 
> files may have something but I was unable to find by myself.
> 
>  
> 
> Thanks in advance. Any help will be appreciated.
> 
>  
> 
> -
> 
> MSc.   Rodrigo Aluizio
> 
> Centro de Estudos do Mar/UFPR
> Laboratório de Micropaleontologia
> Avenida Beira Mar s/n - CEP 83255-000
> Pontal do Paraná - PR - Brasil
> Fone: (0**41) 3455-3620 ramal 217
> Fax: (0**41) 3455-3623
> 
> 
> 
> 
>   [[alternative HTML version deleted]]
> 
> 

The contents of this e-mail are confidential and may be subject to legal 
privilege.
 If you are not the intended recipient you must not use, disseminate, 
distribute or
 reproduce all or any part of this e-mail or attachments.  If you have received 
this
 e-mail in error, please notify the sender and delete all material pertaining 
to this
 e-mail.  Any opinion or views expressed in this e-mail are those of the 
individual
 sender and may not represent those of The New Zealand Institute for Plant and
 Food Research Limited.

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


Re: [R] Plotting with plot()

2009-03-23 Thread Peter Dalgaard

Satindra Chakravorty wrote:

I am trying to use the plot() function in R but can't seem to generate any
plots.  I am running R version 2.7.1 on a server machine to which I am
remotely connected on a Linux platform.

I typed

plot(trees$volume)

as a test and got back the ">" prompt with no plot in my current window or
in any other window.

I read in another post on this list about opening a graphics device with
x11().  I tried that and got an ERROR message:
Error in X11
   unable to start device X11
In x11() : unable to open connection to X11 display ' '.

Please let me know what I should do.


Either look for Rplots.pdf in your current remote directory (after 
leaving R), or make sure that you have an X11 display to which R can 
connect and you can view remotely. The latter depends on how you are 
connected and from which kind of system. If you are using SSH from 
another Linux system running X11, it could be as simple as remembering 
the -Y option to ssh. If you are running Windows, things are more 
complicated (VNC is one option).



--
   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
~~ - (p.dalga...@biostat.ku.dk)  FAX: (+45) 35327907

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


[R] two different date formats in the same variable

2009-03-23 Thread Farrel Buchinsky
How does one convert to a date format when survey respondents have
used two different date formats whilst entering their data. There were
clearly told to use mm/dd/ but humans being humans some entered
mm/dd/yy. There was even validity checks on the forms but I allowed
them to be overridden since the data is more holy than the format.

The data was downloaded as a csv and read.csv was used to read in.
There are several date variables (for example date of birth, date of
diagnosis). Some became character vectors and others become factor
vectors. Nevertheless I have accomplished most of what I want using
lines such as

strptime(init.consent$consent.rec,"%m/%d/%Y")
strptime(x,"%m/%d/%Y")
 as.Date(x, "%m/%d/%Y")

But what happens when a few of the entries get messed up because they
are actually formatted %m/%d/%y.

Is there a robust date formatter? Alternatively how would one code
(presumably using regular expressions) a transforamtion or
substitution only on the errant entries and thereby turn 06/25/04 into
06/25/2004 and 03/03/59 into 03/03/1959?

sessionInfo()
R version 2.8.1 (2008-12-22)
i386-pc-mingw32


Farrel Buchinsky

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


Re: [R] If statement generates two outputs

2009-03-23 Thread Carl Witthoft

>From: Wacek Kusnierczyk 
>Date: Sun, 22 Mar 2009 22:58:49 +0100


>just for fun, you could do this with multiassignment, e.g., using the 
>(highly experimental and premature!) rvalues:


>source('http://miscell.googlecode.com/svn/rvalues/rvalues.r') 
>if (TRUE)


>   c(df1, df2) := list(4:8, 9:13)

>dput(df1)
># 4:8
>dput(df2)
># 9:13


---
Now THAT's what I call an overloaded operator!   ^_^

But seriously:  can someone explain to me what's going on in the 
rvalues.r code?  I tried a simple experiment, replacing ":=" with a 
"colec" in the code, and of course the line


c(df1, df2) colec list(4:8, 9:13)


just gives me a "syntax error" response.   Clearly I need a pointer to 
some documentation about how the colon and equals sign get "special 
treatment somewhere inside R.


thanks
Carl

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


Re: [R] variance/mean

2009-03-23 Thread Wacek Kusnierczyk

(this post suggests a patch to the sources, so i allow myself to divert
it to r-devel)

Bert Gunter wrote:
> x a numeric vector, matrix or data frame. 
> y NULL (default) or a vector, matrix or data frame with compatible
> dimensions to x. The default is equivalent to y = x (but more efficient). 
>
>   
bert points to an interesting fragment of ?var:  it suggests that
computing var(x) is more efficient than computing var(x,x), for any x
valid as input to var.  indeed:

set.seed(0)
x = matrix(rnorm(1), 100, 100)

library(rbenchmark)
benchmark(replications=1000, columns=c('test', 'elapsed'),
   var(x),
   var(x, x))
#test elapsed
# 1var(x)   1.091
# 2 var(x, x)   2.051

that's of course, so to speak, unreasonable:  for what var(x) does is
actually computing the covariance of x and x, which should be the same
as var(x,x). 

the hack is that if y is given, there's an overhead of memory allocation
for *both* x and y when y is given, as seen in src/main/cov.c:720+.
incidentally, it seems that the problem can be solved with a trivial fix
(see the attached patch), so that

set.seed(0)
x = matrix(rnorm(1), 100, 100)

library(rbenchmark)
benchmark(replications=1000, columns=c('test', 'elapsed'),
   var(x),
   var(x, x))
#test elapsed
# 1var(x)   1.121
# 2 var(x, x)   1.107

with the quick checks

all.equal(var(x), var(x, x))
# TRUE
   
all(var(x) == var(x, x))
# TRUE

and for cor it seems to make cor(x,x) slightly faster than cor(x), while
originally it was twice slower:

# original
benchmark(replications=1000, columns=c('test', 'elapsed'),
   cor(x),
   cor(x, x))
#test elapsed
# 1cor(x)   1.196
# 2 cor(x, x)   2.253
   
# patched
benchmark(replications=1000, columns=c('test', 'elapsed'),
   cor(x),
   cor(x, x))
#test elapsed
# 1cor(x)   1.207
# 2 cor(x, x)   1.204

(there is a visible penalty due to an additional pointer test, but it's
10ms on 1000 replications with 1 data points, which i think is
negligible.)

> This is as clear as I would know how to state. 

i believe bert is right.

however, with the above fix, this can now be rewritten as:

"
x: a numeric vector, matrix or data frame. 
y: a vector, matrix or data frame with dimensions compatible to those of x. 
By default, y = x. 
"

which, to my simple mind, is even more clear than what bert would know
how to state, and less likely to cause the sort of confusion that
originated this thread.

the attached patch suggests modifications to src/main/cov.c and
src/library/stats/man/cor.Rd.
it has been prepared and checked as follows:

svn co https://svn.r-project.org/R/trunk trunk
cd trunk
# edited the sources
svn diff > cov.diff
svn revert -R src
patch -p0 < cov.diff

tools/rsync-recommended
./configure
make
make check
bin/R
# subsequent testing within R

if you happen to consider this patch for a commit, please be sure to
examine and test it carefully first.

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


Re: [R] two different date formats in the same variable

2009-03-23 Thread jim holtman
Try using 'strsplit' to split your string on the '/' and then create a
series of 'if's to determine how you want to output the new string.
You will probably need this approach since you may have to check the
validity and ranges of the numbers.

On Mon, Mar 23, 2009 at 8:06 PM, Farrel Buchinsky  wrote:
> How does one convert to a date format when survey respondents have
> used two different date formats whilst entering their data. There were
> clearly told to use mm/dd/ but humans being humans some entered
> mm/dd/yy. There was even validity checks on the forms but I allowed
> them to be overridden since the data is more holy than the format.
>
> The data was downloaded as a csv and read.csv was used to read in.
> There are several date variables (for example date of birth, date of
> diagnosis). Some became character vectors and others become factor
> vectors. Nevertheless I have accomplished most of what I want using
> lines such as
>
> strptime(init.consent$consent.rec,"%m/%d/%Y")
> strptime(x,"%m/%d/%Y")
>  as.Date(x, "%m/%d/%Y")
>
> But what happens when a few of the entries get messed up because they
> are actually formatted %m/%d/%y.
>
> Is there a robust date formatter? Alternatively how would one code
> (presumably using regular expressions) a transforamtion or
> substitution only on the errant entries and thereby turn 06/25/04 into
> 06/25/2004 and 03/03/59 into 03/03/1959?
>
> sessionInfo()
> R version 2.8.1 (2008-12-22)
> i386-pc-mingw32
>
>
> Farrel Buchinsky
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>



-- 
Jim Holtman
Cincinnati, OH
+1 513 646 9390

What is the problem that you are trying to solve?

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


Re: [R] error statement: missing value where TRUE/FALSE needed

2009-03-23 Thread jim holtman
I would assume the expression in the 'if' is giving an NA as a result:

> if (NA) 1 else 2
Error in if (NA) 1 else 2 : missing value where TRUE/FALSE needed
>


On Mon, Mar 23, 2009 at 10:38 AM, Michael Curran  wrote:
> Hi list,
>
> I want to try Gibbs sampling as a method of estimating a markov-switching
> model of a mean-deviating, pth-order autoregressive process with time
> varying transition probabilities via R and am using a code originally
> written by another person; I attach the useful pdf document explaining the
> code. When I run the code, I get an error message:
>
> Error in if (r < vQ[i]) { : missing value where TRUE/FALSE needed
>
> I am using R.2.6.0 on Windows XP. Can anyone tell me what the error message
> means?
>
> The specific line in the code is:
>
> if (r The code is available at: http://www.michael-curran.com/gibbs.html
> and the two data sets are available at:
> http://www.geocities.jp/atsmatsumoto/ci.txt
> and
> http://www.geocities.jp/atsmatsumoto/callrate.txt
>
> Note: I saved the first data set as ci.txt and the second as boj.txt and so
> these are the file names that the code loads. If anyone uses the code and
> the data and manages to get it to work, I would gratefully appreciate if
> they could tell me what adjustments they made to the code in order to do so.
>
> Kind regards,
>
> Michael
> --
> Michael Curran
> Candidate for the MPhil in Economics
> Cambridge University
> http://www.michael-curran.com/
>
>
>
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
>



-- 
Jim Holtman
Cincinnati, OH
+1 513 646 9390

What is the problem that you are trying to solve?

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


[R] Summarizing each row into a frequency table

2009-03-23 Thread Daren Tan
I have a matrix containing -1, 0, 1, however certain rows will not
have all 3 numbers. I have written some codes to compute the frequency
table of how many -1s, 0s, 1s per row, but it is very ugly and not
efficient if there are more than 3 numbers. Please suggest.

m <- rbind(sample(0:1, replace=T, 10), sample(-1:1, replace=T, 10))
m.table <- t(apply(m, 1, function(x) c(sum(x==-1, na.rm=T), sum(x==0,
na.rm=T), sum(x==1, na.rm=T)) ))
m.table <- prop.table(m.table, 1)*100
colnames(m.table) <- -1:1

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


Re: [R] How to set up a function for "Central Limit Theorem"

2009-03-23 Thread j verzani
Greg Snow  imail.org> writes:

> 
> There is also the clt.examp function in the TeachingDemos package.
> 


An interactive demo can also be found at http://www.math.csi.cuny.edu/gWidgetsWWW/ex-clt.html>http://www.math.csi.cuny.edu/gWidgetsWWW/ex-clt.html.

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


  1   2   >