Re: [R] Sampling problem

[wangwallace]

Michael, I really appreciate your help.

but I got the following error message when I wan trying to run the function
written by you:

Error in out[i, ] <- apply(help[, c(grp1 + 1, grp2 + 5)], 2, sample, 1) : 
  number of items to replace is not a multiple of replacement length

I am not quite sure why would this happen.

As a novice of R, these functions are kinda complex for me. I am wondering
if it is doable without using loops like that.

Again, thank you so much!!!  
-- 
View this message in context: 
http://r.789695.n4.nabble.com/Sampling-problem-tp3043804p3044249.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] simulate survival data using median survival time

[Kere Klein]

Dear All,

I like to know how to simulate survival data using median (or mean) survival 
time. Any help will be greatly appreciated.

Best wishes,
Kere

Kerenaftali Klein PhD| Biostatistician | Queensland Clinical Trials & 
Biostatistics Centre 
The University of Queensland | School of Population Health | Building 33, Level 
1| Princess Alexandra Hospital |Ipswich Road | Woolloongabba QLD 4102 | 
Australia Ph: +61 7 3176 3062| Fax: +61 7 3176 6826 | Email: k.kl...@uq.edu.au 
| Web: http://www.sph.uq.edu.au/qctbc
__
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] Sampling problem

[Michael Bedward]

On 16 November 2010 16:10, wangwallace  wrote:
>
> Michael, I really appreciate your help.
>
> but I got the following error message when I wan trying to run the function
> written by you:
>
> Error in out[i, ] <- apply(help[, c(grp1 + 1, grp2 + 5)], 2, sample, 1) :
>  number of items to replace is not a multiple of replacement length

Did the data.frame or matrix you were sampling have the same general
form as the example you posted previously ?  Can you give me a small
example that causes the error ?

> I am not quite sure why would this happen.
>
> As a novice of R, these functions are kinda complex for me. I am wondering
> if it is doable without using loops like that.

I wasn't sure exactly what you wanted so the function was meant to be
general and easy to modify. It is often possible to use constructs
other than loops in R, though that doesn't mean the code will always
be either faster or clearer. But you'll need to describe your
requirements in more precise terms (short, clear examples are good)
for folks here to suggest methods.

>
> Again, thank you so much!!!

No worries. If you can provide an example that generates the error we
should be able to get further.

Michael

__
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] Integrate to 1? (gauss.quad)

[Enrico Schumann]

 
gauss--hermite weights the function f(x) to be integrated by exp(-nodes^2),
ie, it uses the 'trick' that exp(nodes^2) * exp(-nodes^2) = 1. so your sum
should be exp(qq$nodes^2) * exp(-qq$nodes^2) * f(qq$nodes^2) =
exp(qq$nodes^2) * qq$weights * f(qq$nodes). hence

sum(exp(qq$nodes^2) * (1/(s*sqrt(2*pi))) * exp(-((qq$nodes-mu)^2/(2*s^2))) *
qq$weights)

or 

sum(exp(qq$nodes^2) * qq$weights * myNorm(qq$nodes))

should work. 

hth,
enrico

> -Ursprüngliche Nachricht-
> Von: r-help-boun...@r-project.org 
> [mailto:r-help-boun...@r-project.org] Im Auftrag von Doran, Harold
> Gesendet: Montag, 15. November 2010 18:11
> An: Douglas Bates
> Cc: r-help@r-project.org
> Betreff: Re: [R] Integrate to 1? (gauss.quad)
> 
> Thank you, Doug. I am still missing something here. Should 
> this simply be sum(f(x_i) * w_i) where x_i is node i and w_i 
> is the weight at node i? So, my function f(x) = 
> (1/(s*sqrt(2*pi)))  * exp(-((qq$nodes-mu)^2/(2*s^2))) would 
> then be multiplied only by the weights, qq$weights
> 
> Here is some reproducible code:
> 
> > library(statmod)
> > Q <- 5
> > mu <- 0; s <- 1
> > qq <- gauss.quad(Q, kind = 'hermite')
> > sum((1/(s*sqrt(2*pi)))  * exp(-((qq$nodes-mu)^2/(2*s^2))) * 
> > qq$weights)
> [1] 0.5775682
> 
> > sum(qq$weights)
> [1] 1.772454
> 
> > (1/(s*sqrt(2*pi)))  * exp(-((qq$nodes-mu)^2/(2*s^2)))
> [1] 0.05184442 0.25198918 0.39894228 0.25198918 0.05184442
> 
> > dnorm(qq$nodes)
> [1] 0.05184442 0.25198918 0.39894228 0.25198918 0.05184442
> 
> > -Original Message-
> > From: dmba...@gmail.com [mailto:dmba...@gmail.com] On Behalf Of 
> > Douglas Bates
> > Sent: Sunday, November 14, 2010 2:28 PM
> > To: Doran, Harold
> > Cc: r-help@r-project.org
> > Subject: Re: [R] Integrate to 1? (gauss.quad)
> > 
> > I don't know about the statmod package and the gauss.quad 
> function but 
> > generally the definition of Gauss-Hermite quadrature is 
> with respect 
> > to the function that is multiplied by exp(-x^2) in the 
> integrand.  So 
> > your example would reduce to summing the weights.
> > 
> > On Sun, Nov 14, 2010 at 11:18 AM, Doran, Harold 
>  wrote:
> > > Does anyone see why my code does not integrate to 1?
> > >
> > > library(statmod)
> > >
> > > mu <- 0
> > > s <- 1
> > > Q <- 5
> > > qq <- gauss.quad(Q, kind='hermite')
> > > sum((1/(s*sqrt(2*pi)))  * exp(-((qq$nodes-mu)^2/(2*s^2))) * 
> > > qq$weights)
> > >
> > > ### This does what's it is supposed to myNorm <- function(theta) 
> > > (1/(s*sqrt(2*pi)))  * exp(-((theta-mu)^2/(2*s^2))) 
> integrate(myNorm, 
> > > -Inf, Inf) __
> > > 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.
> No virus found in this incoming message.
> Checked by AVG - www.avg.com
> Version: 9.0.869 / Virus Database: 271.1.1/3256 - Release 
> Date: 11/14/10 08:34:00
> 

__
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] simulate survival data using median survival time

[Michael Bedward]

Hi Kere,

Step 1 is choose an appropriate distribution :)  Do you have one in
mind ?  Or are you interested in examining the effects of survival
times having the same mean but generated with alternative
distributions ?

One ready-rolled alternative is the SimSurv method in package prodlim.
To find (many) others you can use the "sos" package...

library(sos)
findFn("simulate survival data")

Which returns 66 matches when I just tried it.

Michael

On 16 November 2010 18:04, Kere Klein  wrote:
> Dear All,
>
> I like to know how to simulate survival data using median (or mean) survival 
> time. Any help will be greatly appreciated.
>
> Best wishes,
> Kere
>
> Kerenaftali Klein PhD| Biostatistician | Queensland Clinical Trials & 
> Biostatistics Centre
> The University of Queensland | School of Population Health | Building 33, 
> Level 1| Princess Alexandra Hospital |Ipswich Road | Woolloongabba QLD 4102 | 
> Australia Ph: +61 7 3176 3062| Fax: +61 7 3176 6826 | Email: 
> k.kl...@uq.edu.au | Web: http://www.sph.uq.edu.au/qctbc
> __
> 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] problem with PDF/postcript, cannot change paper size: "‘mode(width)’ and ‘mod e(height)’ differ between new and previous "

[Emmanuel Levy]

Hi,

The pdf function would not let me change the paper size and gives me
the following warning:

pdf("figure.pdf", width="6", height="10")
Warning message:
‘mode(width)’ and ‘mode(height)’ differ between new and previous
 ==> NOT changing ‘width’ & ‘height’

If I use the option   paper = "a4r",   it does not give me a warning
but still prints on a square region (it does not use the entire page).

Two days ago I updated R and associated packages. I'm not sure if this
could be the cause? Even if I restart a new session I keep getting the
same error.  I also tried X11.options(reset=TRUE) and it does not
help.

Thanks for any hint you could give me, this is really annoying (the
function postscript gives the same error so I cannot make any
figure!).
all the best,

Emmanuel

> sessionInfo()
R version 2.12.0 (2010-10-15)
Platform: x86_64-pc-linux-gnu (64-bit)

locale:
 [1] LC_CTYPE=en_CA.UTF-8   LC_NUMERIC=C
 [3] LC_TIME=en_CA.UTF-8LC_COLLATE=en_CA.UTF-8
 [5] LC_MONETARY=C  LC_MESSAGES=en_CA.UTF-8
 [7] LC_PAPER=en_CA.UTF-8   LC_NAME=C
 [9] LC_ADDRESS=C   LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_CA.UTF-8 LC_IDENTIFICATION=C

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

__
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] glmer, Error: Downdated X'X is not positive definite,49

[Annika]

Dear list,
 
I am new to this list and I am new to the world of R. Additionally I am not
very firm in statistics either but have to deal. So here is my problem:
I have a dataset (which I attach at the end of the post) with a binomial
response variable (alive or not) and three fixed factors
(trapping,treat,sex). I do have repeated measures and would like to include
one (enclosure) random factor.
I chose to use the glmer function (and tried the lmer as well). 
  
m1<-glmer(alive~treat*sex*trapping*ID+(1|enclosure),family=binomial,data=survadult)
Both give the following error message:
   Error in mer_finalize(ans) : Downdated X'X is not positive definite,
49.
   Additionally: Warning:
   glm.fit: fitted probabilities numerically 0 or 1 occurred 
I get that there is something wrong in my dataframe, but I don`t understand
what and how I could change it.
If I run a glm without the random factor
   m1<-glm(alive~treat*sex*trapping*ID,family=binomial,data=survadult)
it works out really fine, but is the wrong analysis in my opinion. (Does glm
actually consider the repeated measures? Does glmer?)

Is glmer the right function for me and how do I deal with the error I get?

Thanks for any help for a non-professional :)

Annika


trappingtreat   sex alive   ID  enclosure
trapping1   m/m f   1   400 1
trapping1   m/m f   1   401 1
trapping1   m/m f   1   402 1
trapping1   m/m m   0   403 1
trapping1   m/m m   1   404 1
trapping1   m/m m   1   405 1
trapping1   m/p f   1   406 2
trapping1   m/p f   1   407 2
trapping1   m/p f   1   408 2
trapping1   m/p m   1   409 2
trapping1   m/p m   1   410 2
trapping1   m/p m   1   411 2
trapping1   p/p f   1   412 3
trapping1   p/p f   1   413 3
trapping1   p/p f   1   414 3
trapping1   p/p m   1   415 3
trapping1   p/p m   1   416 3
trapping1   p/p m   1   417 3
trapping1   p/p f   1   418 4
trapping1   p/p f   1   419 4
trapping1   p/p f   1   420 4
trapping1   p/p m   1   421 4
trapping1   p/p m   1   422 4
trapping1   p/p m   1   423 4
trapping1   m/m f   1   424 5
trapping1   m/m f   1   425 5
trapping1   m/m f   1   426 5
trapping1   m/m m   1   427 5
trapping1   m/m m   1   428 5
trapping1   m/m m   1   429 5
trapping1   p/m f   1   430 6
trapping1   p/m f   1   431 6
trapping1   p/m f   1   432 6
trapping1   p/m m   1   433 6
trapping1   p/m m   1   434 6
trapping1   p/m m   1   435 6
trapping1   m/p f   1   436 7
trapping1   m/p f   1   437 7
trapping1   m/p f   1   438 7
trapping1   m/p m   1   439 7
trapping1   m/p m   1   440 7
trapping1   m/p m   1   441 7
trapping1   p/p f   1   442 8
trapping1   p/p f   1   443 8
trapping1   p/p f   1   444 8
trapping1   p/p m   1   445 8
trapping1   p/p m   1   446 8
trapping1   p/p m   1   447 8
trapping1   m/m f   1   448 9
trapping1   m/m f   1   449 9
trapping1   m/m f   1   450 9
trapping1   m/m m   1   451 9
trapping1   m/m m   1   452 9
trapping1   m/m m   1   453 9
trapping1   p/m f   1   454 10
trapping1   p/m f   1   455 10
trapping1   p/m f   1   456 10
trapping1   p/m m   1   457 10
trapping1   p/m m   1   458 10
trapping1   p/m m   1   459 10
trapping1   m/m f   1   460 11
trapping1   m/m f   1   461 11
trapping1   m/m f   1   462 11
trapping1   

Re: [R] problem with PDF/postcript, cannot change paper size: "‘mode(width)’ and ‘mod e(height)’ differ between new and previous "

[Emmanuel Levy]

Update - sorry for the stupid question, let's say it's pretty late.

For those who may be as tired as I am and get the same warning, the
paper size should be given as an integer!


On 16 November 2010 04:17, Emmanuel Levy  wrote:
> Hi,
>
> The pdf function would not let me change the paper size and gives me
> the following warning:
>
> pdf("figure.pdf", width="6", height="10")
> Warning message:
> ‘mode(width)’ and ‘mode(height)’ differ between new and previous
>         ==> NOT changing ‘width’ & ‘height’
>
> If I use the option   paper = "a4r",   it does not give me a warning
> but still prints on a square region (it does not use the entire page).
>
> Two days ago I updated R and associated packages. I'm not sure if this
> could be the cause? Even if I restart a new session I keep getting the
> same error.  I also tried X11.options(reset=TRUE) and it does not
> help.
>
> Thanks for any hint you could give me, this is really annoying (the
> function postscript gives the same error so I cannot make any
> figure!).
> all the best,
>
> Emmanuel
>
>> sessionInfo()
> R version 2.12.0 (2010-10-15)
> Platform: x86_64-pc-linux-gnu (64-bit)
>
> locale:
>  [1] LC_CTYPE=en_CA.UTF-8       LC_NUMERIC=C
>  [3] LC_TIME=en_CA.UTF-8        LC_COLLATE=en_CA.UTF-8
>  [5] LC_MONETARY=C              LC_MESSAGES=en_CA.UTF-8
>  [7] LC_PAPER=en_CA.UTF-8       LC_NAME=C
>  [9] LC_ADDRESS=C               LC_TELEPHONE=C
> [11] LC_MEASUREMENT=en_CA.UTF-8 LC_IDENTIFICATION=C
>
> attached base packages:
> [1] stats     graphics  grDevices utils     datasets  methods   base
>

__
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] Who uses R with multicore, SNOW or CUDA package for resource intense computing?

[erik]

Hi,

who of you in this forum uses R (http://www.r-project.org/) with the
multicore, SNOW or CUDA packages, so for advanced calculations that need
more power than a workstation CPU? On which hardware do you compute these
scripts? At home/ at work or do you have data center access somewhere?

The background of these questions is the following: I am currently writing
my M.Sc. thesis about R and High-Performance-Computing and need a strong
knowledge about who actually uses R. I read that R had 1 million users in
2008, bu thats more or less the only user statistics I could find on this
topic - so I hope for your answers!

Sincerely Heinrich
-- 
View this message in context: 
http://r.789695.n4.nabble.com/Who-uses-R-with-multicore-SNOW-or-CUDA-package-for-resource-intense-computing-tp3044485p3044485.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] odfWeave - "Format error discovered in the file in sub-document content.xml at 2, 4047 (row, col)"

[Søren Højsgaard]

When using odfWeave on an OpenOffice input document, I can not open the output 
document. I get the message

"Format error discovered in the file in sub-document content.xml at 2,4047 
(row,col)"

Can anyone help me on this? (Apologies if this has been discussed before; I 
have not been able to find any info...)

Info:
I am using R.2.12.0 on Windows 7 (64 bit). I have downloaded the XML package 
from http://www.stats.ox.ac.uk/pub/RWin/bin/windows/contrib/2.12/ and I have 
compiled odfWeave myself

Best regards
Søren

[[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] odfWeave - "Format error discovered in the file in sub-document content.xml at 2, 4047 (row, col)"

[Mike Marchywka]

From: soren.hojsga...@agrsci.dk
To: r-h...@stat.math.ethz.ch
Date: Tue, 16 Nov 2010 11:32:06 +0100
Subject: [R] odfWeave - "Format error discovered in the file in sub-document 
content.xml at 2, 4047 (row, col)"


When using odfWeave on an OpenOffice input document, I can not open the output 
document. I get the message

"Format error discovered in the file in sub-document content.xml at 2,4047 
(row,col)"

Can anyone help me on this? (Apologies if this has been discussed before; I 
have not been able to find any info...)

well, if it really means line 2 you could post the first few lines. Did you 
expect a line
with 4047 columns? 




Info:

e from http://www.stats.ox.ac.uk/pub/RWin/bin/windows/contrib/2.12/ and I have 
compiled odfWeave myself

Best regards
Søren

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

[R] LATTICE. On skip, index.cond with a formula like Y~X|A+B

[ottorino]

Dear invaluable R-list,
my present problem is arranging/removing some panels in a lattice plot.

Please consider the following:

df.data <-
cbind.data.frame(expand.grid(SUBJ=1:5,
 TREAT=LETTERS[1:4],
 REF=letters[1:4]
 )
 )
df.data <-
df.data[order(df.data$TREAT,df.data$REF), ]
df.data <-
cbind.data.frame(df.data[,2:3],
X=1:5,
Y= c(
seq(1,100,length=5),
seq(1,10,length=5),
seq(1,50,length=5),
1:5,
##
seq(1,10,length=5),
seq(1,100,length=5),
seq(1,50,length=5),
1:5,
##
seq(1,10,length=5),
seq(1,50,length=5),
seq(1,100,length=5),
1:5,
##
seq(1,10,length=5),
seq(1,50,length=5),
1:5,
seq(1,100,length=5)
)
)

# Plot the data:
# Note the interaction in the formula
# with this formula lattice draws a strip in the 
# form A:B in one line
(my.plot <-
xyplot(Y ~ X|TREAT:REF,
   data = df.data,
   layout = c(4,4),
   aspect=1,
   panel = function(x, y, ...){
   panel.xyplot(x, y, lty=1, ...)
   panel.abline(coef(lm(y ~ x), col=2))
   })
 )

# Now imagine that I only need the elements off the diagonal,
# say the inferior ones.
# So I can write some lines using a convenient call to 
# "skip" and "index.cond"

update(my.plot,
   skip= c(0,0,0,1,0,0,1,1),
   index.cond=list(c(2,3,4,7,8,12)),
   layout=c(3,3))

# The lines above draw ALMOST what I'm looking for.
# My real problem is that I have factors with long 
# and unshrinkable names, that do not fit easily in one line.
#
# For sake of clarity and for other reasons I would like to plot
# something like the one below, with two strips and with the shingle
# BUT removing all the panels except the off diagonal inferior elements.

(my.plotTwoLines <-
 xyplot(Y ~ X|TREAT+REF,
   data = df.data,
   layout = c(4,4),
   aspect=1,
   strip= strip.custom(stile=3),
   panel = function(x, y, ...){
   panel.xyplot(x, y, lty=1, ...)
   panel.abline(coef(lm(y ~ x), col=2))
   })
 )

# The formula |TREAT+REF results in a 
# "index.cond" which is obviously a list with one element for
# each factor

my.plotTwoLines$index.cond

# In this case I think I can only remove columns or rows, but not single
# panels. 
# Any suggestion ?

# Thanks in advance for your time.

-- 
Ottorino-Luca Pantani, Università di Firenze
Dip.to di Scienze delle Produzioni Vegetali,
del Suolo e dell'Ambiente Forestale (DiPSA)
P.zle Cascine 28 50144 Firenze Italia
Ubuntu 10.04 -- GNOME 2.30.2
GNU Emacs 23.1.50.1 (x86_64-pc-linux-gnu, GTK+ Version
2.18.0)
ESS version 5.11 -- R 2.12.0

__
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] Offset in glm poisson using R vs Exposure in Stata

[Columbine Caroline Waring]

R-helpers,

I am hoping to find someone who uses both R and program Stata for GLMs.

I am a beginner R user, finding my own way through; learning code etc. at the 
same time as learning the statistics I need to complete my project.

What I have is the code from Stata  and am trying to reproduce the same 
analysis in R - my program of choice.

. glm count md ms rf sg, family(poisson)
exposure(effort) eform



I am lost at the point of finding the equivalent code for 'exposure'. 

Having looked at a few forums and 'googled'. I thought 'offset', used as
offset=(log(Eff))   or the equivalent+offset(log(Eff))   would produce the 
desired effect.

Incidentally my code was:
glm(Count~md+ms+rf+sg+offset(Eff),family=poisson,data=DepthHabGen)

(Making use of glm{stats})

However, offset does not seem to be equivalent to 'exposure' in Stata. As 
coefficients and log likelhood estimates differ.

So I asked the following questions:

1. Do both programs produce the same results without 'exposure' i.e. glm models 

Yes,  log likelihoods and coefficients are the same.

2. How about using the unintuitive non logged " offset=Eff" ?

Coefficients and log likelihoods still differ.

So now I defer to those with more experience, does anyone know how to produce 
the equivalent of Stata's glm  'exposure' in R?

Thank you in advance,
Columbine
  
[[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] Debugging segfault in foreach

[Duncan Murdoch]

On 15/11/2010 10:53 PM, Steve Lianoglou wrote:

Hi,

I'm using R-2.12 on a linux 64bit machine.

When I run a chunk of code inside a foreach() %do% { ...} or %dopar%
{...} (with doMC backend) I keep getting a segfault. Running the
*same* code within lapply(something, function(x) ... ) doesn't result
in any segfaults. I'll paste the output below, but I'm not sure it
would be helpful.

I'm more curious how to go about smoking out what's causing this, as
getting back "into" the code block that is crashing is tricky inside
the code block executed by foreach().

If anybody has any suggestions on how they debug with foreach, I'd be
happy to hear them. Setting .veborse=TRUE isn't providing anything
informative.

Anyway, I'll paste the segfault bomb below, followed by sessionInfo()


I have very little experience debugging on Linux, but presumably what 
you want to do is to use a post-mortem debugger to analyze the core dump 
(no idea how to do that), or run R under gdb or another debugger to 
catch things at the time the segfault occurs.  There's a command line 
option to do that; try R --help to see it.  (But it won't give you help 
on gdb, you'll need other sources of info for that.)


Duncan Murdoch



 Segfault dump ==

result of evaluating expression:

  *** caught segfault ***
address 0x30, cause 'memory not mapped'

Traceback:
  1: format(x[[i]], ..., justify = justify)
  2: format.data.frame(x, digits = digits, na.encode = FALSE)
  3: as.matrix(format.data.frame(x, digits = digits, na.encode = FALSE))
  4: print.data.frame(r)
  5: print(r)
  6: doTryCatch(return(expr), name, parentenv, handler)
  7: tryCatchOne(expr, names, parentenv, handlers[[1L]])
  8: tryCatchList(expr, classes, parentenv, handlers)
  9: tryCatch({repeat {args<- nextElem(it)if
(obj$verbose) {cat(sprintf("evaluation # %d:\n", i))
  print(args)}for (a in names(args)) assign(a\
, args[[a]], pos = envir, inherits = FALSE)r<-
tryCatch(eval(expr, envir = envir), error = function(e) e)if
(obj$verbose) {cat("result of evaluating expression\
:\n")print(r)}
tryCatch(accumulator(list(r), i), error = function(e) {
cat("error calling combine function:\n")print(e)
  NULL}) \
i<- i + 1}}, error = function(e) {if
(!identical(conditionMessage(e), "StopIteration"))
stop(simpleError(conditionMessage(e), expr))})
10: doSEQ(obj, substitute(ex), parent.frame())
11: foreach(chr = chrs, .packages = "GenomicFeaturesX", .verbose =
TRUE) %do% {.gc<- duplicate(gcache, pre.load = NULL)
   on.exit(dispose(.gc))cat("===", chr, "===\n")   \
  apaSummary(expt, .gc, chr, pvd.policy, utr.index, polya.variants,
 gene.collapse)}
12: apaSummaryCrank(bpm[[1]], gcr, chrs[21:22], gene.collapse = "longest")

Possible actions:
1: abort (with core dump, if enabled)
2: normal R exit
3: exit R without saving workspace
4: exit R saving workspace

= sessionInfo ==

R version 2.12.0 (2010-10-15)
Platform: x86_64-unknown-linux-gnu (64-bit)

locale:
  [1] LC_CTYPE=en_US.UTF-8   LC_NUMERIC=C
  [3] LC_TIME=en_US.UTF-8LC_COLLATE=en_US.UTF-8
  [5] LC_MONETARY=C  LC_MESSAGES=en_US.UTF-8
  [7] LC_PAPER=en_US.UTF-8   LC_NAME=C
  [9] LC_ADDRESS=C   LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C

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

other attached packages:
  [1] reshape2_1.0  ShortRead_1.8.0   lattice_0.19-13
  [4] ggplot2_0.8.8 proto_0.3-8   reshape_0.8.3
  [7] doMC_1.2.1multicore_0.1-3   foreach_1.3.0
[10] codetools_0.2-2   iterators_1.0.3   data.table_1.5.1
[13] plyr_1.2.1GenomeGraphs_1.10.0   biomaRt_2.6.0
[16] bitops_1.0-4.1Rsamtools_1.2.0   RSQLite_0.9-2
[19] DBI_0.2-5 Biostrings_2.18.0 GenomicFeaturesX_0.2
[22] GenomicFeatures_1.2.0 GenomicRanges_1.2.1   IRanges_1.8.2

loaded via a namespace (and not attached):
  [1] annotate_1.28.0  AnnotationDbi_1.12.0 Biobase_2.10.0
  [4] BSgenome_1.18.0  hwriter_1.2  RCurl_1.4-3
  [7] rtracklayer_1.10.2   stringr_0.4  XML_3.2-0
[10] xtable_1.5-6

-steve



__
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] scatterplot with filled circles

[DrCJones]

for a simple scatterplot: 

plot(X ~ Y, type = 'p', col = 'red') 

this produces red-edged circles, but I want to fill in the circles.

can this be done? I checked '? plot' already but couldn't find what I was
looking for. 

cheers

 



-- 
View this message in context: 
http://r.789695.n4.nabble.com/scatterplot-with-filled-circles-tp3044690p3044690.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] Force evaluation of variable when calling partialPlot

[Tim Howard]

Greg, 
Two thoughts: 
1. It might be possible that 'vars' is a reserved word of sorts and if you 
change the name of your vector RF might be happier
2. A way that works for me is to call importance as follows:
 
sel.imp <- importance(sel.rf, class=NULL, scale=TRUE, type=NULL)
 
and then use the 'names' of the imp data frame to be absolutely clear to RF you 
are talking about the same variables
 
for(i in length(sel.imp){

partialPlot(sel.rf,xdata,names(sel.imp[i]),which.class=1,xlab=vars[i],main="")
}
 
 
Hope that helps. 
Tim Howard
 
>>>
Date: Mon, 15 Nov 2010 12:29:08 -0800 (PST)
From: gdillon 
To: r-help@r-project.org
Subject: Re: [R] Force evaluation of variable when calling partialPlot
Message-ID: <1289852948670-3043750.p...@n4.nabble.com>
Content-Type: text/plain; charset=us-ascii


RE: the folloing original question:
> I'm using the randomForest package and would like to generate partial
> dependence plots, one after another, for a variety of variables:
> 
> m <- randomForest( s, ... )
> varnames <- c( "var1", "var2", "var3", "var4" )   # var1..4 are all in
> data frame s
> for( v in varnames ) {
>partialPlot( x=m, pred.data=s, x.var=v )
> }
> 
> ...but this doesn't work, with partialPlot complaining that it can't
> find the variable "v".

I'm having a very similar problem using the partialPlot function. A
simplified version of my code looks like this:

data.in <- paste(basedir,"/sw_climate_dataframe_",root.name,".csv",sep="")
data <- read.csv(data.in)
vars <- c("sm1","precip.spring","tmax.fall","precip.fall")  #selected
variables in data frame
xdata <- as.data.frame(data[,vars])
ydata <- data[,5]
ntree <- 2000

rf.pdplots <- function() {
  sel.rf <- randomForest(xdata,ydata,ntree=ntree,keep.forest=TRUE)
  par(family="sans",mfrow=c(2,2),mar=c(4,3,1,2),oma=c(0,3,0,0),mgp=c(2,1,0))
  for (i in 1:length(vars)) {
print((vars)[i])
partialPlot(sel.rf,xdata,vars[i],which.class=1,xlab=vars[i],main="")
mtext("(Logit of probability of high severity)/2", side=2, line=1,
outer=T)
  }
}

rf.pdplots()

When I run this code, with partialPlots embedded in a function, I get the
following error:
Error in eval(expr, envir, enclos) : object "i" not found

If I just take the code inside the function and run it (not embedded in the
function), it runs just fine. Is there some reason why partialPlots doesn't
like to be called from inside a function?

Other things I've tried/noticed:
1. If I comment out the line with the partialPlots call (and the next line
with mtext), the function runs as expected and prints the variable names one
at a time.
2. If the variable i is defined as a number (e.g., 4) in the global
environment, then the function will run, the names print out one at a time,
and four plots are created. HOWEVER, the plots are all for the last (4th)
variable, BUT the x labels actually are different on each plot (i.e., the
xlab is actually looping through the four values in vars).

Can anyone help me make sense of this? Thanks.

-Greg
-- 

[[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] AOV/LME

[Tal Galili]

Hi janesense,

You might find it useful to go through the tutorial listed here for help:
http://www.r-statistics.com/2010/04/repeated-measures-anova-with-r-tutorials/

Cheers,
Tal

Contact
Details:---
Contact me: tal.gal...@gmail.com |  972-52-7275845
Read me: www.talgalili.com (Hebrew) | www.biostatistics.co.il (Hebrew) |
www.r-statistics.com (English)
--




On Tue, Nov 16, 2010 at 3:31 AM, janesense  wrote:

>
> Hi everyone,
>
> I'm having a little trouble with working out what formula is better to use
> for a repeated measures two way anova. I have two factors, L (five levels)
> and T (two levels). L and T are both crossed factors (all participants do
> all combinations). So, I do:
>
>summary(aov(dat~L*T+Error(participant/(L*T)),data=dat4))
>
> But get different results with:
>
>   anova(lme(dat~L*T,random=~1|participant,data=dat4)).
>
> Rather, the lme results are the same as those with the formula
>
>   summary(aov(dat~L*T+Error(participant),data=dat4))
>
> I know most people advocate the use of lme rather than aov, but I think I
> need to include factors L and T in the error term and can't figure out how
> to do this with lme, so possibly the first aov is more correct. I have a
> balanced design so I don't think it would be a problem to use aov.
>
> The first aov shows an expected significant interaction (unlike lme and the
> second aov), but I don't want to use this if it's incorrect!
>
> If anyone could help, that would be great.
>
>
>
>
> --
> View this message in context:
> http://r.789695.n4.nabble.com/AOV-LME-tp3044118p3044118.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] Debugging segfault in foreach

[Rainer M Krug]

-BEGIN PGP SIGNED MESSAGE-
Hash: SHA1

On 11/16/2010 04:53 AM, Steve Lianoglou wrote:
> Hi,
> 
> I'm using R-2.12 on a linux 64bit machine.
> 
> When I run a chunk of code inside a foreach() %do% { ...} or %dopar%
> {...} (with doMC backend) I keep getting a segfault. Running the
> *same* code within lapply(something, function(x) ... ) doesn't result
> in any segfaults. I'll paste the output below, but I'm not sure it
> would be helpful.

Have you tried using mclapply instead? Does that work? If yes, the
problem is in foreach and I would contact revolution analytics - if no,
the problem is in multicore (or your code ...)

I don't think you can go much further then that, unless you want to
debug foreach (or multicore) to that effect.

Cheers,

Rainer


> 
> I'm more curious how to go about smoking out what's causing this, as
> getting back "into" the code block that is crashing is tricky inside
> the code block executed by foreach().
> 
> If anybody has any suggestions on how they debug with foreach, I'd be
> happy to hear them. Setting .veborse=TRUE isn't providing anything
> informative.
> 
> Anyway, I'll paste the segfault bomb below, followed by sessionInfo()
> 
>  Segfault dump ==
> 
> result of evaluating expression:
> 
>  *** caught segfault ***
> address 0x30, cause 'memory not mapped'
> 
> Traceback:
>  1: format(x[[i]], ..., justify = justify)
>  2: format.data.frame(x, digits = digits, na.encode = FALSE)
>  3: as.matrix(format.data.frame(x, digits = digits, na.encode = FALSE))
>  4: print.data.frame(r)
>  5: print(r)
>  6: doTryCatch(return(expr), name, parentenv, handler)
>  7: tryCatchOne(expr, names, parentenv, handlers[[1L]])
>  8: tryCatchList(expr, classes, parentenv, handlers)
>  9: tryCatch({repeat {args <- nextElem(it)if
> (obj$verbose) {cat(sprintf("evaluation # %d:\n", i))
>  print(args)}for (a in names(args)) assign(a\
> , args[[a]], pos = envir, inherits = FALSE)r <-
> tryCatch(eval(expr, envir = envir), error = function(e) e)if
> (obj$verbose) {cat("result of evaluating expression\
> :\n")print(r)}
> tryCatch(accumulator(list(r), i), error = function(e) {
> cat("error calling combine function:\n")print(e)
>  NULL}) \
>i <- i + 1}}, error = function(e) {if
> (!identical(conditionMessage(e), "StopIteration"))
> stop(simpleError(conditionMessage(e), expr))})
> 10: doSEQ(obj, substitute(ex), parent.frame())
> 11: foreach(chr = chrs, .packages = "GenomicFeaturesX", .verbose =
> TRUE) %do% {.gc <- duplicate(gcache, pre.load = NULL)
>   on.exit(dispose(.gc))cat("===", chr, "===\n")   \
>  apaSummary(expt, .gc, chr, pvd.policy, utr.index, polya.variants,
> gene.collapse)}
> 12: apaSummaryCrank(bpm[[1]], gcr, chrs[21:22], gene.collapse = "longest")
> 
> Possible actions:
> 1: abort (with core dump, if enabled)
> 2: normal R exit
> 3: exit R without saving workspace
> 4: exit R saving workspace
> 
> = sessionInfo ==
> 
> R version 2.12.0 (2010-10-15)
> Platform: x86_64-unknown-linux-gnu (64-bit)
> 
> locale:
>  [1] LC_CTYPE=en_US.UTF-8   LC_NUMERIC=C
>  [3] LC_TIME=en_US.UTF-8LC_COLLATE=en_US.UTF-8
>  [5] LC_MONETARY=C  LC_MESSAGES=en_US.UTF-8
>  [7] LC_PAPER=en_US.UTF-8   LC_NAME=C
>  [9] LC_ADDRESS=C   LC_TELEPHONE=C
> [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
> 
> attached base packages:
> [1] grid  stats graphics  grDevices utils datasets  methods
> [8] base
> 
> other attached packages:
>  [1] reshape2_1.0  ShortRead_1.8.0   lattice_0.19-13
>  [4] ggplot2_0.8.8 proto_0.3-8   reshape_0.8.3
>  [7] doMC_1.2.1multicore_0.1-3   foreach_1.3.0
> [10] codetools_0.2-2   iterators_1.0.3   data.table_1.5.1
> [13] plyr_1.2.1GenomeGraphs_1.10.0   biomaRt_2.6.0
> [16] bitops_1.0-4.1Rsamtools_1.2.0   RSQLite_0.9-2
> [19] DBI_0.2-5 Biostrings_2.18.0 GenomicFeaturesX_0.2
> [22] GenomicFeatures_1.2.0 GenomicRanges_1.2.1   IRanges_1.8.2
> 
> loaded via a namespace (and not attached):
>  [1] annotate_1.28.0  AnnotationDbi_1.12.0 Biobase_2.10.0
>  [4] BSgenome_1.18.0  hwriter_1.2  RCurl_1.4-3
>  [7] rtracklayer_1.10.2   stringr_0.4  XML_3.2-0
> [10] xtable_1.5-6
> 
> -steve
> 


- -- 
Rainer M. Krug, PhD (Conservation Ecology, SUN), MSc (Conservation
Biology, UCT), Dipl. Phys. (Germany)

Centre of Excellence for Invasion Biology
Natural Sciences Building
Office Suite 2039
Stellenbosch University
Main Campus, Merriman Avenue
Stellenbosch
South Africa

Tel:+33 - (0)9 53 10 27 44
Cell:   +27 - (0)8 39 47 90 42
Fax (SA):   +27 - (0)8 65 16 27 82
Fax (D) :   +49 - (0)3 21 21 25 22 44
Fax (FR):   +33 - (0)9 58 10 27 44
email:  rai...@krugs.de

Skype:  RMkrug
-BEGIN PGP SIGNATURE-
Vers

Re: [R] scatterplot with filled circles

[B.-Markus Schuller]

Hey,

the parameter in the plot command you are looking for is "pch".
For filled circles type:

> plot(X ~ Y, type = 'p', col = 'red', pch=16)

The 16 determines the filled circles.
Cheers,

Mango

-
B.-Markus Schuller aka Mango

Sensory Ecology Group
Max-Planck-Institute for Ornithology
82319 Seewiesen, Germany

phone: +49 (0)8157 932 -378
fax:   +49 (0)8157 932 -344
email: schul...@orn.mpg.de
http://www.orn.mpg.de/nwg/abtsiemers.html
-
Never run for the bus.
Never skip tea.


On 16.11.2010 13:28, DrCJones wrote:


for a simple scatterplot:

plot(X ~ Y, type = 'p', col = 'red')

this produces red-edged circles, but I want to fill in the circles.

can this be done? I checked '? plot' already but couldn't find what I was
looking for.

cheers







__
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] scatterplot with filled circles

[Ben Bolker]

DrCJones  gmail.com> writes:

> for a simple scatterplot: 
> 
> plot(X ~ Y, type = 'p', col = 'red') 
> 
> this produces red-edged circles, but I want to fill in the circles.
 
plot(X ~ Y, type = 'p', col = 'red', pch=16) 

  Are you sure you want X~Y and not Y~X?

__
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] Offset in glm poisson using R vs Exposure in Stata

[Ben Bolker]

Columbine Caroline Waring  hotmail.com> writes:

> I am hoping to find someone who uses both R and program Stata for GLMs.
[snip]
> What I have is the code from Stata  and am trying to reproduce the same
analysis in R - my program of choice.
> 
> . glm count md ms rf sg, family(poisson)
> exposure(effort) eform
> 
> I am lost at the point of finding the equivalent code for 'exposure'. 
> 
> Having looked at a few forums and 'googled'. I thought 'offset', used as   
offset=(log(Eff))   or the
> equivalent+offset(log(Eff))   would produce the desired effect.
> 
> Incidentally my code was:  
  glm(Count~md+ms+rf+sg+offset(Eff),family=poisson,data=DepthHabGen)
> 
> (Making use of glm{stats})
> 

  Based on your discussion above did you mean

glm(count~md+ms+rf+sg+offset(log(Eff)),
   family=poisson,data=DepthHabGen)

 ? is this just a typo ?

according to http://www.stata.com/help.cgi?glm

 exposure(varname) include ln(varname) in model with coefficient
  constrained to 1
 offset(varname)   include varname in model with coefficient
  constrained to 1


> However, offset does not seem to be equivalent to 'exposure' in Stata. As
coefficients and log likelhood
> estimates differ.
> 
> So I asked the following questions:
> 
> 1. Do both programs produce the same results without
> 'exposure' i.e. glm models 
> 
> Yes,  log likelihoods and coefficients are the same.
> 
> 2. How about using the unintuitive non logged " offset=Eff" ?
> 
> Coefficients and log likelihoods still differ.

  You're not doing anything obviously wrong.  Are you sure your
"effort" and "Eff" variables are the same, i.e. nothing got mangled
moving to R?

  I don't use Stata, perhaps someone else can try.  Posting a small
reproducible example would be helpful.

__
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] scatterplot with filled circles

[DrCJones]

Fantastic thanks!
-- 
View this message in context: 
http://r.789695.n4.nabble.com/scatterplot-with-filled-circles-tp3044690p3044772.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] scatterplot with filled circles

[DrCJones]

"Are you sure you want X~Y and not Y~X? "

You are right of course - sorry that was just a typo ;) 
-- 
View this message in context: 
http://r.789695.n4.nabble.com/scatterplot-with-filled-circles-tp3044690p3044776.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] odfWeave - "Format error discovered in the file in sub-document content.xml at 2, 4047 (row, col)"

[Søren Højsgaard]

Dear Mike,

Good point - thanks. The lines that caused the error mentioned above are simply:

<<>>=
x <- 1:10
x
@

I could add that the document 'simple.odt' (which comes with odfWeave) causes 
the same error - but at row=109, col=1577

> sessionInfo()
R version 2.12.0 (2010-10-15)
Platform: x86_64-pc-mingw32/x64 (64-bit)

locale:
[1] LC_COLLATE=Danish_Denmark.1252  LC_CTYPE=Danish_Denmark.1252
LC_MONETARY=Danish_Denmark.1252 LC_NUMERIC=C
LC_TIME=Danish_Denmark.1252

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

other attached packages:
[1] MASS_7.3-8  odfWeave_0.7.14 XML_3.2-0.1 lattice_0.19-13

loaded via a namespace (and not attached):
[1] tools_2.12.0


Regards
Søren

-Oprindelig meddelelse-
Fra: Mike Marchywka [mailto:marchy...@hotmail.com] 
Sendt: 16. november 2010 12:56
Til: Søren Højsgaard; r-h...@stat.math.ethz.ch
Emne: RE: [R] odfWeave - "Format error discovered in the file in sub-document 
content.xml at 2, 4047 (row, col)"









From: soren.hojsga...@agrsci.dk
To: r-h...@stat.math.ethz.ch
Date: Tue, 16 Nov 2010 11:32:06 +0100
Subject: [R] odfWeave - "Format error discovered in the file in sub-document 
content.xml at 2, 4047 (row, col)"


When using odfWeave on an OpenOffice input document, I can not open the output 
document. I get the message

"Format error discovered in the file in sub-document content.xml at 2,4047 
(row,col)"

Can anyone help me on this? (Apologies if this has been discussed before; I 
have not been able to find any info...)

well, if it really means line 2 you could post the first few lines. Did you 
expect a line
with 4047 columns? 




Info:
I am using R.2.12.0 on Windows 7 (64 bit). I have downloaded the XML package 
from http://www.stats.ox.ac.uk/pub/RWin/bin/windows/contrib/2.12/ and I have 
compiled odfWeave myself

Best regards
Søren

[[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] scatterplot with filled circles

[Dennis Murphy]

Try adding  pch = 16 to your plot call.

HTH,
Dennis

On Tue, Nov 16, 2010 at 4:28 AM, DrCJones wrote:

>
> for a simple scatterplot:
>
> plot(X ~ Y, type = 'p', col = 'red')
>
> this produces red-edged circles, but I want to fill in the circles.
>
> can this be done? I checked '? plot' already but couldn't find what I was
> looking for.
>
> cheers
>
>
>
>
>
> --
> View this message in context:
> http://r.789695.n4.nabble.com/scatterplot-with-filled-circles-tp3044690p3044690.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.

[R] Question about readLines

[romzero]

Can i use "readLines" to extract only the linees with a specific word within?

If yes, how?

Tnx for help.
-- 
View this message in context: 
http://r.789695.n4.nabble.com/Question-about-readLines-tp3044701p3044701.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] renaming a list of data frames y calculating with lapply

[David Mendez]

Hello everyone I'm  a new R user and have some questions about lists
indexing, names and the lapply function.

first have a look to my data (names changed) it's a rainfall record for 67
stations, 

str(data)
List of 67
 $ e.somename.xx.txt  :'data.frame':   456 obs. of  5
variables:
  ..$ Y: int [1:456] 1960 1960 1960 1960 1960 1960 1960 1960 1960 1960
...
  ..$ M: int [1:456] 1 2 3 4 5 6 7 8 9 10 ...
  ..$ PRE  : num [1:456] 0 0 0 0 10.6 ...
  ..$ COD  : Factor w/ 2 levels "FE","O": 2 2 2 2 2 2 2 2 2 2 ...
  ..$ piece: num [1:456] -0.874 -0.874 -0.874 -0.874 -0.66 ...
...

on a single a data frame it is possible to use names or <- to change column
name "piece" for something more meaningful  to me like "my_calc"  i don't
know how generalize that to a list of dataframes perhaps something like

lapply(data, names,...)

"piece" is in fact the result of few lines of ugly code with rapply, sapply
and llply the names were lost somewhere...

1. How do i change "piece" name?

I need to subset every dataframe in the list by "$M" and then over this
subset aggregate using mean of "$piece" on a single data frame i do

z <- subset(e999.somedata.xx.txt, M==c(1,2,3))
z <- aggregate(z[,5], list(z[,1]),mean)

I tried this but it doesn't behave the way i expected, (because it creates a
single data frame and a want a list of data frames). 

for  (i in data){
z <- subset(i, i$M==c(1,2,3))
z <- aggregate(z[,5], list(z[,1]),mean)
}

2. how to subset and aggregate?, perhaps lapply again?

a rainfall of thanks to the charitable soul who can help with this basics
questions, i just don get indexation or naming conventions!
-- 
View this message in context: 
http://r.789695.n4.nabble.com/renaming-a-list-of-data-frames-y-calculating-with-lapply-tp3044908p3044908.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] Question about readLines

[David Winsemius]

On Nov 16, 2010, at 7:37 AM, romzero wrote:



Can i use "readLines" to extract only the linees with a specific  
word within?


If yes, how?


Just use:

lines <- readLines(  )
desired <- lines[grep("word", lines]




Tnx for help.
--
View this message in context: 
http://r.789695.n4.nabble.com/Question-about-readLines-tp3044701p3044701.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.


David Winsemius, MD
West Hartford, CT

__
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] Question about readLines

[jim holtman]

x <- readLines('yourFile')
x <- x[grepl('myWord', x)]


On Tue, Nov 16, 2010 at 7:37 AM, romzero  wrote:
>
> Can i use "readLines" to extract only the linees with a specific word within?
>
> If yes, how?
>
> Tnx for help.
> --
> View this message in context: 
> http://r.789695.n4.nabble.com/Question-about-readLines-tp3044701p3044701.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.
>



-- 
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] Breslow-Day test

[Robert Ruser]

Dear R Users,
I'm looking for a package that allows to test hypothesis about a
homogeneity of odds ratio in k 2x2 tables. I know that Breslow-Day is
suitable but does anybody could me point out a package? I found diffR,
but as far as I see this package is for IRT theory.

Best,
Robert

__
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 article in forbes

[Mark Leeds]

Here's an article in Forbes about R in case anyone's interested.

http://blogs.forbes.com/smcnally/2010/11/10/names-you-need-to-know-in-2011-r-data-analysis-software/

[[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] Breslow-Day test

[David Winsemius]

On Nov 16, 2010, at 10:08 AM, Robert Ruser wrote:


Dear R Users,
I'm looking for a package that allows to test hypothesis about a
homogeneity of odds ratio in k 2x2 tables. I know that Breslow-Day is
suitable but does anybody could me point out a package? I found diffR,
but as far as I see this package is for IRT theory.


Try searching:
http://search.r-project.org/cgi-bin/namazu.cgi?query=&max=100&result=normal&sort=score&idxname=Rhelp10&idxname=Rhelp08&idxname=Rhelp02&idxname=functions

At least two of the links point to this:

http://www.math.montana.edu/~jimrc/classes/stat524/Rcode/breslowday.test.r

--
David Winsemius, MD
West Hartford, CT

__
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] Re : interpretation of coefficients in survreg AND obtaining the hazard function for an individual given a set of predictors

[Vincent Vinh-Hung]

Thanks for sharing the questions and responses!

Is it possible to appreciate how much the coefficients matter in one
or the other model?
Say, using Biau's example, using coxph, as.factor(grade2 ==
"high")TRUE gives hazard ratio 1.27 (rounded).
As clinician I can grasp this HR as 27% relative increase. I can
relate with other published results.
With survreg the Weibull model gives a coefficient -0.4035245: is it
feasible or meaningful to translate it to HR?

Thanks in advance,

Vincent Vinh-Hung
Radiation Oncology,
Geneva University Hospitals

On Sun, Nov 14, 2010 at 6:51 AM, Biau David  wrote:
> Dear R help list,
>
> I am modeling some survival data with coxph and survreg (dist='weibull') using
> package survival. I have 2 problems:
>
> 1) I do not understand how to interpret the regression coefficients in the
> survreg output and it is not clear, for me, from ?survreg.objects how to.
>
> Here is an example of the codes that points out my problem:
> - data is stc1
> - the factor is dichotomous with 'low' and 'high' categories
>
> slr <- Surv(stc1$ti_lr, stc1$ev_lr==1)
>
> mca <- coxph(slr~as.factor(grade2=='high'), data=stc1)
> mcb <- coxph(slr~as.factor(grade2), data=stc1)
> mwa <- survreg(slr~as.factor(grade2=='high'), data=stc1, dist='weibull',
> scale=0)
> mwb <- survreg(slr~as.factor(grade2), data=stc1, dist='weibull', scale=0)
>
>> summary(mca)$coef
> coef
> exp(coef)  se(coef) z  Pr(>|z|)
> as.factor(grade2 == "high")TRUE 0.2416562  1.273356 0.2456232
> 0.9838494  0.3251896
>
>> summary(mcb)$coef
>   coef exp(coef)
> se(coef) z Pr(>|z|)
> as.factor(grade2)low -0.2416562 0.7853261 0.2456232 -0.9838494
> 0.3251896
>
>> summary(mwa)$coef
> (Intercept) as.factor(grade2 == "high")TRUE
> 7.9068380   -0.4035245
>
>> summary(mwb)$coef
> (Intercept) as.factor(grade2)low
> 7.5033135   0.4035245
>
>
> No problem with the interpretation of the coefs in the cox model. However, i 
> do
> not understand why
> a) the coefficients in the survreg model are the opposite (negative when the
> other is positive) of what I have in the cox model? are these not the log(HR)
> given the categories of these variable?

No. survreg() fits accelerated failure models, not proportional
hazards models.   The coefficients are logarithms of ratios of
survival times, so a positive coefficient means longer survival.


> b) how come the intercept coefficient changes (the scale parameter does not
> change)?

Because you have reversed the order of the factor levels.  The
coefficient of that variable changes sign and the intercept changes to
compensate.


> 2) My second question relates to the first.
> a) given a model from survreg, say mwa above, how should i do to extract the
> base hazard and the hazard of each patient given a set of predictors? With the
> hazard function for the ith individual in the study given by  h_i(t) =
> exp(\beta'x_i)*\lambda*\gamma*t^{\gamma-1}, it doesn't look like to me that
> predict(mwa, type='linear') is \beta'x_i.

No, it's beta'x_i for the accelerated failure parametrization of the
Weibull.  In terms of the CDF

F_i(t) = F_0( exp((t+beta'x_i)/scale) )

So you need to multiply by the scale parameter and change sign to get
the log hazard ratios.


> b) since I need the coefficient intercept from the model to obtain the scale
> parameter  to obtain the base hazard function as defined in Collett
> (h_0(t)=\lambda*\gamma*t^{\gamma-1}), I am concerned that this coefficient
> intercept changes depending on the reference level of the factor entered in 
> the
> model. The change is very important when I have more than one predictor in the
> model.

As Terry Therneau pointed out recently in the context of the Cox
model, there is no such thing as "the" baseline hazard.  The baseline
hazard is the hazard when all your covariates are equal to zero, and
this depends on how you parametrize.  In mwa, zero is grade2="low", in
mwb, zero is grade2="high", so the hazard at zero has to be different
in the two cases.

 -thomas

--
Thomas Lumley
Professor of Biostatistics
University of Auckland

__
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.dendrogram() plot margins

[Uwe Ligges]

On 15.11.2010 18:35, francois fauteux wrote:

Hello,

Is it possible to remove those extra margins on the "sample" axis from
plot.dendrogram:

par(oma=c(0,0,0,0),mar=c(0,0,0,0))
ddr<-as.dendrogram(hclust(dist(matrix(sample(1:1000,200),nrow=100
stats:::plot.dendrogram(ddr,horiz=F,axes=F,yaxs="i",leaflab="none")



I guess you want

 par(xaxs="i")
 stats:::plot.dendrogram(ddr,horiz=F,axes=F,yaxs="i",leaflab="none")

Uwe Ligges





vs.

stats:::plot.dendrogram(ddr,horiz=T,axes=F,yaxs="i",leaflab="none")

What variable / line of code corresponds to this additional margin space? I
would like to modify the code to remove the extra space and have that margin
equal to that when horiz=T, for plotting multiple dendrograms for one images
on the same device.

Thanks in advance, best.

[[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] Breslow-Day test

[Marc Schwartz]

On Nov 16, 2010, at 9:18 AM, David Winsemius wrote:

> 
> On Nov 16, 2010, at 10:08 AM, Robert Ruser wrote:
> 
>> Dear R Users,
>> I'm looking for a package that allows to test hypothesis about a
>> homogeneity of odds ratio in k 2x2 tables. I know that Breslow-Day is
>> suitable but does anybody could me point out a package? I found diffR,
>> but as far as I see this package is for IRT theory.
> 
> Try searching:
> http://search.r-project.org/cgi-bin/namazu.cgi?query=&max=100&result=normal&sort=score&idxname=Rhelp10&idxname=Rhelp08&idxname=Rhelp02&idxname=functions
> 
> At least two of the links point to this:
> 
> http://www.math.montana.edu/~jimrc/classes/stat524/Rcode/breslowday.test.r


You might also want to look at:

  http://statweb.calpoly.edu/aschaffn/418/documents/Script8.01.pdf

which has code that includes the Tarone adjustment based upon:

On Heterogeneity Tests Based on Efficient Scores
Tarone RE 
Biometrika 72, 91-95. 1985

Supported by Breslow in:
Statistics in Epidemiology: The Case-Control Study
N. E. Breslow
Journal of the American Statistical Association
Vol. 91, No. 433 (Mar., 1996), pp. 14-28

and also look at the code for the woolf() function in ?mantelhaen.test

HTH,

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

Re: [R] [klaR package] [NaiveBayes] warning message numerical 0 probability

[Fabon Dzogang]

Thank you very much for this clarification,

Sincerely,
Fabon Dzogang.

2010/11/15 Uwe Ligges :
>
>
> On 03.11.2010 16:26, Fabon Dzogang wrote:
>>
>> Hi,
>>
>> I run R 2.10.1 under ubuntu 10.04 LTS (Lucid Lynx) and klaR version 0.6-4.
>>
>> I compute a model over a 2 classes dataset (composed of 700 examples).
>> To that aim, I use the function NaiveBayes provided in the package
>> klaR.
>> When I then use the prediction function : predict(my_model, new_data).
>> I get the following warning :
>>
>> "In FUN(1:747[[747L]], ...) : Numerical 0 probability with observation
>> 458"
>>
>> As I did not find any documentation or any discussion concerning this
>> warning message, I looked in the klaR source code and found the
>> following line in predict.NaiveBayes.R :
>>
>> "warning("Numerical 0 probability with observation ", i)"
>
>
> Within Naive Bayes in order to calculate the posteriori probabilities of
> the classes it is necessary to calculate the probabilities of the
> observations given the classes. The function NaiveBayes prints a warning
> if all these probabilities are numerical 0, i.e. that the observation
> has a numerical probability of 0 for *all* classes. Usually this is only
> the case when the obs. is an extreme outlier.
>
> I will change the warning to say "all classes" in further releases of klaR.
>
> Best wishes,
> Uwe Ligges
>
>
>
>> Unfortunately, it is hard to get a clear picture of the whole process
>> reading the code. I wonder if someone could help me with the meaning
>> of this warning message.
>>
>> Sorry I did not provide an example, but I could not simulate the same
>> message over a small toy example.
>>
>> Thank you,
>>
>> Fabon Dzogang.
>>
>> __
>> 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.
>



-- 
Fabon Dzogang

__
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] Simple error handling in R

[Aleksey Naumov]

Hi R experts,

I am looking for a simple error handling approach, whereby I could stop
function execution with a customized error message. For example:

for (i in 1:10) {
   if (i == 5)
  # I'd like to be able to stop right here with an error message I have
complete control over
}

The problem with stop() is that I cannot control completely what gets
printed to the terminal, even with stop(call.=FALSE) there is the "Error:"
string. I've worked through examples (whatever few there are) for try() and
tryCatch() and I still cannot understand how to do this. If I supply my own
error handler function with tryCatch(..., error=function(e) ...) I can
control the error message, but the loop continues on to i=6, etc.

So I am struggling with error handling in R... It' seems its a lot simpler
and more consistent e.g. in Python.

Any help would be greatly appreciated!

Aleksey

[[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] Integrating functions / vector arithmetic

[Eduardo de Oliveira Horta]

Thanks, I'll try that.

Regarding your question "Why are you using assignments to indicate the
return values of functions?",

I'd probably answer that it is a matlab vice... But I'll change my habits
soon enough!

Best regards,

Eduardo Horta

On Tue, Nov 16, 2010 at 12:27 PM, Erich Neuwirth <
erich.neuwi...@univie.ac.at> wrote:

> Some remarks:
> Why are you using assignments to indicate the return values of functions?
> This is R idiom:
> > f<-function(u){
> >  f<-0
> >  for (j in 1:4){
> >f<-f+coeff[j]*(zeta(j)(u))
> >  }
> >  f
> > }
> >
>
> If you only want the inner product of 2 vector, "outer" probably is an
> overkill.
>
> g <- function(u){
>  sum(coeff * zeta(1:4)(u))
> }
>
> should be enough.
>
> And plot(g) does not work because g is not vectorized.
>
> plot(Vectorize(g))
>
> will work
>
> On Nov 16, 2010, at 1:18 AM, Eduardo de Oliveira Horta wrote:
>
> > Hello,
> >
> > I was trying to build some functions which I would like to integrate over
> an
> > interval using the function 'integrate' from the 'stats' package. As an
> > example, please consider the function
> >
> > h(u)=sin(pi*u) + sqrt(2)*sin(pi*2*u) + sqrt(3)*sin(pi*3*u) +
> 2*sin(pi*4*u)
> >
> > Two alternative ways to 'build' this function are as in f and g below:
> >
> > coeff<-sqrt(1:4)
> >
> > zeta<-function(i) {force(i); function(u){
> >  zeta<-sqrt(2)+sin(pi*i*u)
> > }}
> >
> > f<-function(u){
> >  f<-0
> >  for (j in 1:4){
> >f<-f+coeff[j]*(zeta(j)(u))
> >  }
> >  f<-f
> > }
> >
> > g<-function(u){
> >  g<-crossprod(coeff,zeta(1:4)(u))
> > }
> >
> > Obviously, f and g are equivalent, but in the actual code I am writing, g
> is
> > much faster. However, I can't seem to integrate (neither to plot) g. In
> > fact, these are the error messages I get:
> >
> >> print(f(.1))
> > [1] 4.443642
> >> print(g(.1))
> > [,1]
> > [1,] 4.443642
> >
> > Everything ok until here... but...
> >
> > When using plot(), I get this:
> >
> >> plot(f) #plots, as expected.
> >> plot(g)
> > Error in crossprod(coeff,zeta(1:4)(u)) : arguments not compatible
> > Besides that: Warning message:
> > In pi * i * u :
> >  longer object length is not a multiple of shorter object length
> >
> > When using integrate(), I get this:
> >
> >> integrate(f,0,1)
> > 1.004172 with absolute error < 2.5e-13
> >> integrate(g,0,1)
> > Error in crossprod(coeff,zeta(1:4)(u)) : arguments not compatible
> > Besides that: Warning message:
> > In pi * i * u :
> >  longer object length is not a multiple of shorter object length
> >>
> >
> > I have already tried some simple 'solutions', for example to set
> > g<-function(u){
> >  vec1<-drop(coeff)
> >  vec2<-drop(zeta(1:4)(u))
> >  g<-drop(crossprod(coeff,zeta(1:4)(u)))
> > }
> >
> > as well as using the %*% operation, but these won't solve my problem.
> >
> > Any suggestions would be welcome. Thanks in advance,
> >
> > Eduardo Horta
> >
> >   [[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.
> >
>
> --
> Erich Neuwirth, University of Vienna
> Faculty of Computer Science
> Center for Computer Science Didactics and Learning Research
> Visit our SunSITE at http://sunsite.univie.ac.at
> Phone: +43-1-4277-39902 Fax: +43-1-4277-39459
>
>

[[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] R article in forbes

[Bogaso Christofer]

Thanks Mark, for sharing such a great article. Really feeling proud that I
am also part of that great community. BTW can anyone share the code on how
to draw the picture displayed there (blur surface plot)?

Thanks,

-Original Message-
From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On
Behalf Of Mark Leeds
Sent: 16 November 2010 20:46
To: r-help@r-project.org
Subject: [R] R article in forbes

Here's an article in Forbes about R in case anyone's interested.

http://blogs.forbes.com/smcnally/2010/11/10/names-you-need-to-know-in-2011-r
-data-analysis-software/

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

[R] Odp: Sampling problem

[Petr PIKAL]

Hi

Here is one way (If I understood what you did ask).

test<-read.table("clipboard", header=T)
> test
  SubID CSE1 CSE2 CSE3 CSE4 WSE1 WSE2 WSE3 WSE4
1 165626224
2 264726623
3 355555545
4 454344452
5 556756441
6 654364373
7 736636521
8 836636547

fff<-function(dat, col1=2, col2=1) {
# col1 are number of columns from fist set and col2 from the second set
sel1<-sample(1:4, col1)
sel2<-sample((1:4)[-sel1], col2)
dat[,c(sel1+1,sel2+5)]
# i presume that your data are same as you posted, if not you has to 
change above values
}

fff(test)
  CSE2 CSE1 WSE3
1562

8634
> fff(test)
  CSE1 CSE2 WSE3
1652

8364
> fff(test)
  CSE1 CSE3 WSE4
1664

8367

If you want to do it 1000 times just use simple loop

result <- vector("list", 1000)
for (i in 1:1000) result[[i]] <- fff(test)

Regards
Petr

r-help-boun...@r-project.org napsal dne 15.11.2010 21:59:21:

> wangwallace  
> Odeslal: r-help-boun...@r-project.org
> 
> 15.11.2010 21:59
> 
> Komu
> 
> r-help@r-project.org
> 
> Kopie
> 
> Předmět
> 
> [R] Sampling problem
> 
> 
> Hey,
> 
> I am hoping someone can help me with a sampling question.
> 
> I have a data frame of 8 variables (the first column is the subjects' 
id):
> 
> SubIDCSE1 CSE2 CSE3 CSE4 WSE1 WSE2 WSE3 WSE4 
>   1  6  5   6   2  6  22   4
>   2  6  4   7   2  6  62   3
>   3  5  5   5   5  5  54   5
>   4  5  4   3   4  4  45   2
>   5  5  6   7   5  6  44   1
>   6  5  4   3   6  4  37   3
>   7  3  6   6   3  6  52   1
>   8  3  6   6   3  6  54   7 

> 
> the 6 variables are categorized into two groups with CSE1, CSE2, CSE3, 
and
> CSE4 in one group and the rest in another group. 
> 
> >sample(data[,2:4],2,replace=FALSE)
> 
>CSE1 CSE2 
> 1  65 
> 2  64 
> 3  55 
> 4  54 
> 5  56 
> 6  54 
> 7  36 
> 8  36 
> 
> Now I want to sample 1 column from another group of variables (i.e., 
WSE1,
> WSE2, WSE3, WSE4), but I want to restrict a vector I am going to sample 
from
> to only those columns that are not correspond to GROUP 1 variables I 
have
> sampled. That is, I want to sample a column from WSE3, WSE4  Columns 
> corresponding to CSE1 and CSE2 (i.e., WSE1, WSE2) need to be dropped. 
> 
> How can I do this? what if I want to repeat this whole process (drawing 
2
> random columns from CSE1, CSE2, CSE3, and CSE4 first, AND then another
> random column from WSE1, WSE2, WSE3, and WSE4) for 1000 times. any 
ideas?
> 
> Many thanks in advance!!
> 
> -- 
> View this message in context: 
http://r.789695.n4.nabble.com/Sampling-problem-
> tp3043804p3043804.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] Simple error handling in R

[David Winsemius]

On Nov 16, 2010, at 10:58 AM, Aleksey Naumov wrote:


Hi R experts,

I am looking for a simple error handling approach, whereby I could  
stop

function execution with a customized error message. For example:

for (i in 1:10) {
  if (i == 5){

#Then insert:
 cat("Informative error message") ; break() }

 # I'd like to be able to stop right here with an error message  
I have

complete control over
}



I think "stop right here" is ambiguous. Are you at a terminal session  
or inside a script or function? Do you want to call a debugger?


And if you were inside a function then just cat("message"); return()


The problem with stop() is that I cannot control completely what gets
printed to the terminal, even with stop(call.=FALSE) there is the  
"Error:"
string. I've worked through examples (whatever few there are) for  
try() and
tryCatch() and I still cannot understand how to do this. If I supply  
my own

error handler function with tryCatch(..., error=function(e) ...) I can
control the error message, but the loop continues on to i=6, etc.

So I am struggling with error handling in R... It' seems its a lot  
simpler

and more consistent e.g. in Python.

Any help would be greatly appreciated!

Aleksey

[[alternative HTML version deleted]]


David Winsemius, MD
West Hartford, CT

__
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] randomForest parameters for image classification

[Deschamps, Benjamin]

I have modified my code since asking my original question. The
classifier is now generated correctly (with a good, low error rate, as
expected). However, I am running into two issues: 

1) I am getting an error at the prediction stage, I get only NA's when I
try to run data down the forest;
2) I run out of memory when generating the forest with more than 200
trees due to the large block of memory already occupied by the training
data

Here is my code:


library(raster)
library(randomForest)

# Set some user variables
fn = "image.pix"
outraster = "output.pix"
training_band = 2
validation_band = 1

# Get the training data
myraster = stack(fn)
training_class = subset(myraster, training_band)
training_class[training_class == 0] = NA
training_class = Which(training_class != 0, cells=TRUE)
training_data = extract(myraster, training_class)
training_response = as.factor(as.vector(training_data[,training_band]))
training_predictors = training_data[,3:nlayers(myraster)]
remove(training_data)

# Create and save the forest
r_tree = randomForest(training_predictors, y=training_response, ntree =
200, keep.forest=TRUE) # Runs out of memory with ntree > ~200
remove(training_predictors, training_response)

# Classify the whole image
predictor_data = subset(myraster, 3:nlayers(myraster))
layerNames(predictor_data) = layerNames(myraster)[3:nlayers(myraster)]
predictions = predict(predictor_data, r_tree, filename=outraster,
format="PCIDSK", overwrite=TRUE, progress="text", type="response") #All
NA!?
remove(predictor_data)


See also a thread I started on
http://stackoverflow.com/questions/4186507/rgdal-efficiently-reading-lar
ge-multiband-rasters about improving the efficiency of collecting the
training data...

Thanks, Benjamin


-Original Message-
From: Liaw, Andy [mailto:andy_l...@merck.com] 
Sent: November 11, 2010 7:02 AM
To: Deschamps, Benjamin; r-help@r-project.org
Subject: RE: [R] randomForest parameters for image classification

Please show us the code you used to run randomForest, the output, as
well as what you get with other algorithms (on the same random subset
for comparison).  I have yet to see a dataset where randomForest does
_far_ worse than other methods.

Andy 

> -Original Message-
> From: r-help-boun...@r-project.org 
> [mailto:r-help-boun...@r-project.org] On Behalf Of Deschamps, Benjamin
> Sent: Tuesday, November 09, 2010 10:52 AM
> To: r-help@r-project.org
> Subject: [R] randomForest parameters for image classification
> 
> I am implementing an image classification algorithm using the
> randomForest package. The training data consists of 31000+ training
> cases over 26 variables, plus one factor predictor variable (the
> training class). The main issue I am encountering is very low overall
> classification accuracy (a lot of confusion between classes). 
> However, I
> know from other classifications (including a regular decision tree
> classifier) that the training and validation data is sound and capable
> of producing good accuracies). 
> 
>  
> 
> Currently, I am using the default parameters (500 trees, mtry not set
> (default), nodesize = 1, replace=TRUE). Does anyone have experience
> using this with large datasets? Currently I need to randomly sample my
> training data because giving it the full 31000+ cases returns 
> an out of
> memory error; the same thing happens with large numbers of 
> trees.  From
> what I read in the documentation, perhaps I do not have 
> enough trees to
> fully capture the training data?
> 
>  
> 
> Any suggestions or ideas will be greatly appreciated.
> 
>  
> 
> Benjamin
> 
> 
>   [[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.
> 
Notice:  This e-mail message, together with any attachme...{{dropped:12}}

__
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] Simple error handling in R

[ottorino]

Il giorno mar, 16/11/2010 alle 10.58 -0500, Aleksey Naumov ha scritto:
> for (i in 1:10) {
>if (i == 5)
>   # I'd like to be able to stop right here with an error message I
> have
> complete control over
> } 

I'm far from expert, but perhaps this is near what you are looking for

for(i in 1:10){
if (i == 5) {browser()}
 print("anything")}

Browser stop the execution and allows you to indagate about

-- 
Ottorino-Luca Pantani, Università di Firenze
Dip.to di Scienze delle Produzioni Vegetali,
del Suolo e dell'Ambiente Forestale (DiPSA)
P.zle Cascine 28 50144 Firenze Italia
Ubuntu 10.04 -- GNOME 2.30.2
GNU Emacs 23.1.50.1 (x86_64-pc-linux-gnu, GTK+ Version
2.18.0)
ESS version 5.11 -- R 2.12.0

__
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 article in forbes/ somewhat OT

[Mike Marchywka]

> From: bogaso.christo...@gmail.com
> To: marklee...@gmail.com; r-help@r-project.org
> Date: Tue, 16 Nov 2010 21:48:10 +0530
> Subject: Re: [R] R article in forbes
>
> Thanks Mark, for sharing such a great article. Really feeling proud that I
> am also part of that great community. BTW can anyone share the code on how
> to draw the picture displayed there (blur surface plot)?

In many cases, the problem is not the R code but actually
getting the data you need to make real-life decisions or form opinions.
I sent a link to one of my rants to the OP on this topic 
but it is a bit political so I will spare the whole list.
I would mention however that many day-to-day uses of R
can be limited due to lack of good data in a computer
readable format- you see people posting here about scraping
stock quotes from human readable web pages but sometimes you need to do that 
with 
data from government sources or they just cater to commercial
products. If you use R at home, you may want to take note
of requests from data sources on ideas about how they can optimize their
data delivery for your needs. 

>
> Thanks,
>
> -Original Message-
> From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On
> Behalf Of Mark Leeds
> Sent: 16 November 2010 20:46
> To: r-help@r-project.org
> Subject: [R] R article in forbes
>
> Here's an article in Forbes about R in case anyone's interested.
>
> http://blogs.forbes.com/smcnally/2010/11/10/names-you-need-to-know-in-2011-r
> -data-analysis-software/
>
> [[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.
  
__
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] rotate column names in large matrix

[Lara Poplarski]

Thank you both, these are very helpful hints.

Chris, could you please suggest how to modify what you sent to also show the
same labels as (horizontal) row names? I have not yet mastered the details
of R graphics...

Many thanks in advance,
Lara

On Mon, Nov 15, 2010 at 3:25 PM, Chris Stubben  wrote:

>
> You could display the matrix as an image then add column names (rotated
> used
> srt) if you really want them.
>
> x <- cor(matrix(rnorm(600), 60, 100))
>
> # set margins with extra space at top and xpd=TRUE to write outside plot
> region
> op<-par(mar=c(1,1,5,1), xpd=TRUE)
>
> # display image without the 90 degree counter clockwise rotation
> image(t(x[nrow(x):1,]), axes=FALSE)
>
> ## add 100 column names
> y<-paste("column", 1:100)
> text( seq(0,1,length=100) , 1.01, y,  pos = 2, srt = 270, offset=0, cex=.7)
>
>
> Chris Stubben
>
>
>
>
>
> Lara Poplarski wrote:
> >
> > I have a large (1600*1600) matrix generated with symnum, that I am using
> > to
> > eyeball the structure of a dataset.
> >
> > I have abbreviated the column names with the abbr.colnames option. One
> way
> > to get an even more compact view of the matrix would be to display the
> > column names rotated by 90 degrees.
> >
> > Any pointers on how to do this would be most useful. Any other tips for
> > displaying the matrix in compact form are of course also welcome.
> >
> >
> >
>
> --
> View this message in context:
> http://r.789695.n4.nabble.com/rotate-column-names-in-large-matrix-tp3043493p3043927.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.

[R] Simple error handling in R

[Aleksey Naumov]

David,

Thank you for your help. I agree, "stop right here" is ambiguous, let me
clarify.

I am at a terminal session, calling my top-level function f(). At  any point
inside f() I'd like to be able to generate an error message, print it to the
terminal and stop execution of f() (but still remain in the R session). No,
I do not need to be able to debug the problem, the error message is meant to
be sufficient for an end user (I am writing it for my office colleagues, and
I don not expect them to do any debugging, merely bring the error message to
my attention).

break would not be sufficient, since it deals with loops, any my use of
for() was just an example illustrating that I'd like to stop function
execution at any point.return() is a good suggestion, except that I would
need to call it from my main function f(), e.g.

f = function(...) {
   # An error condition 1
   cat('an error message')
   return()

   # other error conditions
   # ...
}

My hope would be to wrap error messaging and execution control (i.e.
stopping f()) in one helper function -- to make code more readable and less
verbose, e.g.:

f = function(...) {
   # An error condition 1
   my_error_function(message)

   # ...
}

where my_error_function() simply prints the detailed error message and then
stops execution -- i.e. gets out of f().

Does this make sense?

Thank you,
Aleksey






On Tue, Nov 16, 2010 at 11:16 AM, David Winsemius wrote:

>
> On Nov 16, 2010, at 10:58 AM, Aleksey Naumov wrote:
>
>  Hi R experts,
>>
>> I am looking for a simple error handling approach, whereby I could stop
>> function execution with a customized error message. For example:
>>
>> for (i in 1:10) {
>>  if (i == 5){
>>
> #Then insert:
>  cat("Informative error message") ; break() }
>
>
>  # I'd like to be able to stop right here with an error message I have
>> complete control over
>> }
>>
>>
> I think "stop right here" is ambiguous. Are you at a terminal session or
> inside a script or function? Do you want to call a debugger?
>
> And if you were inside a function then just cat("message"); return()
>
>  The problem with stop() is that I cannot control completely what gets
>> printed to the terminal, even with stop(call.=FALSE) there is the "Error:"
>> string. I've worked through examples (whatever few there are) for try()
>> and
>> tryCatch() and I still cannot understand how to do this. If I supply my
>> own
>> error handler function with tryCatch(..., error=function(e) ...) I can
>> control the error message, but the loop continues on to i=6, etc.
>>
>> So I am struggling with error handling in R... It' seems its a lot simpler
>> and more consistent e.g. in Python.
>>
>> Any help would be greatly appreciated!
>>
>> Aleksey
>>
>>[[alternative HTML version deleted]]
>>
>
> David Winsemius, MD
> West Hartford, CT
>
>

[[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 using function merge_all()

[arjun]

Hi,
   I want to merge 4 data frames with one column in common but I am
getting error message while using this function. Can any one help me out. 
> merge_all(Br,Ki,Lu,Pr,by="Genes")
Error: could not find function "merge_all"

I have installed the package: reshape but I still get this error 
-- 
View this message in context: 
http://r.789695.n4.nabble.com/help-using-function-merge-all-tp3045173p3045173.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] Pass character vector to function argument

[Kevin Ummel]

A bit embarrassed to post this seemingly trivial question, but I can't find 
anything in the archive that's quite relevant:

a1=1
a2=2

obs=objects(pattern=glob2rx("a?"))

I want to utilize 'obs' as a function argument to produce something like:

sum(a1,a2)

Obviously, sum(obs) doesn't work, but I've tried variations of 'eval', 'parse', 
'mget', and 'noquote' without success. What am I missing?

Thanks,
Kevin



[[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] Sampling problem

[wangwallace]

yes, the data.frame is exactly the same as the one I posted earlier. 

I was trying to see if the loop function works. And I got that message. here
below is the syntax I was trying to run, followed by the error message at
the end:
> sampleX<-function(X,nGrp1,nsamples){if(nGrp1>=4)stop("can't sample all
> group 1 variables")
+ out<-matrix(0,nsamples,nGrp1+1)
+ for(i in 1:nsamples){
+ grp1<-sample(4,nGrp1)
+ grp2<-sample((1:4)[-grp1],)
+ out[i,]<-apply(X[,c(grp1+1,grp2+5)],2,sample,1)
+ }
+ out}
> sampleX(help,1,10)
Error in out[i, ] <- apply(X[, c(grp1 + 1, grp2 + 5)], 2, sample, 1) : 
  number of items to replace is not a multiple of replacement length

By the way, it is only a small piece of my data set, which has 12 variables
(or columns) for each group (grp1: CSE1, CSE2, CSE3, CSE4, CSE5, CSE6, CSE7,
CSE8, CSE9, CSE10, CSE11, CSE12; grp2: WSE1, WSE2, WSE3, WSE4, WSE5, WSE6,
WSE7, WSE8, WSE9, WSE10, WSE11, WSE12). I will draw 1000 random samples for
each of the 11 different combinations below: 

combination 1: 1 variable from grp1 + 11 variables from grp2 = 12 variables
combination 2: 2 variable from grp1 + 10 variables from grp2 = 12 variables
combination 3: 3 variable from grp1 + 9 variables from grp2 = 12 variables
combination 4: 4 variable from grp1 + 8 variables from grp2 = 12 variables
combination 5: 5 variable from grp1 + 7 variables from grp2 = 12 variables
combination 6: 6 variable from grp1 + 6 variables from grp2 = 12 variables
combination 7: 7 variable from grp1 + 5 variables from grp2 = 12 variables
combination 8: 8 variable from grp1 + 4 variables from grp2 = 12 variables
combination 9: 9 variable from grp1 + 3 variables from grp2 = 12 variables
combination 10: 10 variable from grp1 + 2 variables from grp2 = 12 variables
combination 11: 11 variable from grp1 + 1 variables from grp2 = 12 variables

As shown above, the sum of the variables in each combination will have to be
12. Also, I want to restrict a vector I am going to sample from to only
those columns that are not correspond to grp1 variables I have sampled. For
example, if I sampled 1 variable, say CSE1, from grp1, the other 11
variables from grp2 should not include WSE1; if I sampled 2 variables, say
CSE1 and CSE2, from grp1, the other 10 variables from grp2 should not
include WSE1 and WSE2.

Anyway, this is a lot more complicated example than the one I described in
my first post. But I think I can modify your function if I wanna apply it to
the large data set with 12 variables for each group, since they basically
share the same method. Now I am wondering where the error message is from.

Again, thanks!! :)


-- 
View this message in context: 
http://r.789695.n4.nabble.com/Sampling-problem-tp3043804p3045095.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] how to do linear regression when dimension is high

[poko2000]

Hi I am a newbie in R.
I have data with dim of 20.
How to use lm if i want to do regression with the whole design matrix? My y
is the first column.
the left are xs.
Thanks a lot.

-- 
View this message in context: 
http://r.789695.n4.nabble.com/how-to-do-linear-regression-when-dimension-is-high-tp3045167p3045167.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] Simple error handling in R

[Aleksey Naumov]

Ottorino,

Thank you for your suggestion. I am not really looking to do any debugging,
but merely to issue an informative error message and stop execution of my
top-level function, called at an R terminal. See me response to David for
more details...

Thank you,
Aleksey

On Tue, Nov 16, 2010 at 11:17 AM, ottorino
wrote:

> Il giorno mar, 16/11/2010 alle 10.58 -0500, Aleksey Naumov ha scritto:
> > for (i in 1:10) {
> >if (i == 5)
> >   # I'd like to be able to stop right here with an error message I
> > have
> > complete control over
> > }
>
> I'm far from expert, but perhaps this is near what you are looking for
>
> for(i in 1:10){
> if (i == 5) {browser()}
>  print("anything")}
>
> Browser stop the execution and allows you to indagate about
>
> --
> Ottorino-Luca Pantani, Università di Firenze
> Dip.to di Scienze delle Produzioni Vegetali,
> del Suolo e dell'Ambiente Forestale (DiPSA)
> P.zle Cascine 28 50144 Firenze Italia
> Ubuntu 10.04 -- GNOME 2.30.2
> GNU Emacs 23.1.50.1 (x86_64-pc-linux-gnu, GTK+ Version
> 2.18.0)
> ESS version 5.11 -- R 2.12.0
>
>

[[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] Pass character vector to function argument

[David Winsemius]

On Nov 16, 2010, at 10:34 AM, Kevin Ummel wrote:

A bit embarrassed to post this seemingly trivial question, but I  
can't find anything in the archive that's quite relevant:


a1=1
a2=2

obs=objects(pattern=glob2rx("a?"))

I want to utilize 'obs' as a function argument to produce something  
like:


sum(a1,a2)



?get
> sum(sapply(obs, get))
[1] 3


Obviously, sum(obs) doesn't work, but I've tried variations of  
'eval', 'parse', 'mget', and 'noquote' without success. What am I  
missing?


Thanks,
Kevin


--
David Winsemius, MD
West Hartford, CT

__
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 vs print ??

[skan]

Hello

What's the differente betwen using "plot" and using "print"  in order to
plot a graph?
For example in order to plot the result of a histogram.

cheers
-- 
View this message in context: 
http://r.789695.n4.nabble.com/plot-vs-print-tp3045256p3045256.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] R article in forbes

[John Fox]

Dear David and Bogasso,

The image is from
, as you'd see
by clicking on the graph in the Forbes article, and it was indeed produced
with scatter3d() in the car package.

Best,
 John


John Fox
Senator William McMaster 
  Professor of Social Statistics
Department of Sociology
McMaster University
Hamilton, Ontario, Canada
web: socserv.mcmaster.ca/jfox


> -Original Message-
> From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org]
On
> Behalf Of David Winsemius
> Sent: November-16-10 11:30 AM
> To: Bogaso Christofer
> Cc: r-help@r-project.org
> Subject: Re: [R] R article in forbes
> 
> 
> On Nov 16, 2010, at 11:18 AM, Bogaso Christofer wrote:
> 
> > Thanks Mark, for sharing such a great article. Really feeling proud
> > that I am also part of that great community. BTW can anyone share the
> > code on how to draw the picture displayed there (blur surface plot)?
> >
> 
> Again... a small bit of searching at Baron's site (only 2 hits for s
search
> "prestige eduction income 3d":
> 
> http://finzi.psych.upenn.edu/R/library/car/html/scatter3d.html
> 
> Turns out none of them is exactly that plot, but this is pretty close
(after
> click-drag-rotating the interactive RGL display)
> 
> require(car)
>  if(interactive() && require(rgl) && require(mgcv))
{scatter3d(prestige ~
> income + education, fit=c( "additive"),
>   data=Prestige)
>   }
> 
> (screengrab of X11 window attached.)

__
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] Pass character vector to function argument

[Ivan Calandra]

Hi,

It would have been nice to see what you've tried.

Here is what I got (maybe not the easiest, but looks like it works)
a1 <- 1
a2 <- 2
obs <- objects(pattern=glob2rx("a?"))
sum(unlist(mget(obs, envir=.GlobalEnv)))
[1] 3

It also works with vectors:
a1 <- 1:10
a2 <- 11:20
sum(unlist(mget(obs, envir=.GlobalEnv)))
[1] 210

Is it what you're looking for?

HTH,
Ivan



Le 11/16/2010 16:34, Kevin Ummel a écrit :

A bit embarrassed to post this seemingly trivial question, but I can't find 
anything in the archive that's quite relevant:

a1=1
a2=2

obs=objects(pattern=glob2rx("a?"))

I want to utilize 'obs' as a function argument to produce something like:

sum(a1,a2)

Obviously, sum(obs) doesn't work, but I've tried variations of 'eval', 'parse', 
'mget', and 'noquote' without success. What am I missing?

Thanks,
Kevin



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



--
Ivan CALANDRA
PhD Student
University of Hamburg
Biozentrum Grindel und Zoologisches Museum
Abt. Säugetiere
Martin-Luther-King-Platz 3
D-20146 Hamburg, GERMANY
+49(0)40 42838 6231
ivan.calan...@uni-hamburg.de

**
http://www.for771.uni-bonn.de
http://webapp5.rrz.uni-hamburg.de/mammals/eng/mitarbeiter.php

__
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] Revolutions Blog: October Roundup

[David Smith]

I write about R every weekday at the Revolutions blog:
 http://blog.revolutionanalytics.com
and every month I post a summary of articles from the previous month
of particular interest to readers of r-help.

In case you missed them, here are some articles related to R from the
month of October:

Reviews of the winners and finalists of the 2010 ggplot2 case study
competition: http://bit.ly/ckJxHZ

We have published a new article "R is Hot", with interviews from a
dozen R users in industry and academia: http://bit.ly/d7rZI9

A new code highlighting tool for displaying R code on the web:
http://bit.ly/bra4ap

A white paper and video describing how an anomaly in the 2000 US
Census data is revealed using the RevoScaleR package:
http://bit.ly/dzFXgk

A suggestion for a workflow for R projects to promote transparency,
maintainability, modularity, portability, reproducibility and
efficiency: http://bit.ly/bopvbR

A video of using the "code snippets" feature in the Revolution R IDE
to insert templates of R code: http://bit.ly/dgJ7aR

A new competition to build a predictive model to identify popular R
packages in CRAN: http://bit.ly/cL9NI7

A profile of Rhipe (Hadoop/R integration) author and new Revolution
employee Saptarshi Guha: http://bit.ly/9k7ABg .

Saptarshi Guha's presentation at Hadoop World on using Rhipe to
analyze VOIP quality data was profiled in the SD Times:
http://bit.ly/a81qdy

A lattice chart illustrates the impact of the new Google Instant on
paid search: http://bit.ly/cSi7d8

Improvements in R 2.12.0, released on October 15: http://bit.ly/9nti8L

Three upcoming R courses from Statistics.com: http://bit.ly/a1Xbdo

R is nominated for (http://bit.ly/cZUl99) and wins
(http://bit.ly/9PdD5C) a major Open Source award in New Zealand.

An article in InfoWorld notes that R is "a programming language on the
rise": http://bit.ly/aLLdU9

O'Reilly has published a "Rough Cuts" preview of the forthcoming "R
Cookbook" by Paul Teetor: http://bit.ly/cUC5da

Revolution Analytics names Lee Edlefsen as Chief Scientist: http://bit.ly/a9b7zc

Other non-R-related stories in the past month included the Data
Science Venn Diagram (http://bit.ly/d5AzgN), busting gay stereotypes
with data analysis (http://bit.ly/aRza7H), World Statistics Day
(http://bit.ly/agxfjO), Arthur C Clarke's uncanny predictions from
1964, an article in the NYT about the language of Statistics
(http://bit.ly/cMpCXH), SAS's battle against open source
(http://bit.ly/abP8sz), and a Tufte map of the US economic stimulus.

On a lighter note, we had: why we should get rid of pennies
(http://bit.ly/avS0TJ), the Möbius Bagel (http://bit.ly/clPLCm), and a
groan-worthy Physics pun (http://bit.ly/bzGOjY).

There are new R user groups in Toronto (http://bit.ly/bWhJyw), Houston
(http://bit.ly/c0XFGp) and Cincinnati/Dayton (http://bit.ly/97FpZx).

The R Community Calendar has also been updated at:
http://blog.revolutionanalytics.com/calendar.html

If you're looking for more articles about R, you can find summaries
from previous months at http://blog.revolutionanalytics.com/roundups/.
Join the Revolution mailing list at
http://revolutionanalytics.com/newsletter to be alerted to new
articles on a monthly basis.

As always, thanks for the comments and please keep sending suggestions
to me at da...@revolutionanalytics.com . Don't forget you can also
follow the blog using an RSS reader like Google Reader, or by
following me on Twitter (I'm @revodavid).

# David Smith

--
David M Smith 
VP of Marketing, Revolution Analytics  http://blog.revolutionanalytics.com
Tel: +1 (650) 330-0553 x205 (Palo Alto, CA, USA)

__
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 do linear regression when dimension is high

[Douglas Bates]

On Tue, Nov 16, 2010 at 10:30 AM, poko2000  wrote:

> Hi I am a newbie in R.
> I have data with dim of 20.
> How to use lm if i want to do regression with the whole design matrix? My y
> is the first column.
> the left are xs.
> Thanks a lot.

Do you have the data stored in a matrix or in a data frame?  If it is
in a data frame and the first column is named "y" then you can use a
call like

lm(y ~ ., mydata)

The '.' on the right hand side of the formula is expanded to a formula
incorporating of the columns in the data frame except for the
response.  For example

> summary(lm(Fertility ~ ., swiss))

Call:
lm(formula = Fertility ~ ., data = swiss)

Residuals:
 Min   1Q   Median   3Q  Max
-15.2743  -5.2617   0.5032   4.1198  15.3213

Coefficients:
 Estimate Std. Error t value Pr(>|t|)
(Intercept)  66.91518   10.70604   6.250 1.91e-07
Agriculture  -0.172110.07030  -2.448  0.01873
Examination  -0.258010.25388  -1.016  0.31546
Education-0.870940.18303  -4.758 2.43e-05
Catholic  0.104120.03526   2.953  0.00519
Infant.Mortality  1.077050.38172   2.822  0.00734

Residual standard error: 7.165 on 41 degrees of freedom
Multiple R-squared: 0.7067, Adjusted R-squared: 0.671
F-statistic: 19.76 on 5 and 41 DF,  p-value: 5.594e-10

__
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] Force evaluation of variable when calling partialPlot

[Greg Dillon]

Tim,
Thanks for the suggestions. I tried them both, and while neither solved my 
problem directly, they helped me to get to a solution.

What I've realized is that partialPlots, for some reason, always looks to 
the Global environment when evaluating the "x.var" argument. So any 
variables created within a function (such as "i" in my looping example, or 
even "sel.imp" when I tried your 2nd suggestion) are not seen. From what I 
can tell, this only seems to be true for the "x.var" argument. Other 
arguments seem to be able to take variables from within the function's 
environment/namespace. My solution was to just modify my original function 
to kick a copy of "i" out to Global, as follows:

rf.pdplots <- function() {
  sel.rf <- randomForest(xdata,ydata,ntree=ntree,keep.forest=TRUE)
 
par(family="sans",mfrow=c(2,2),mar=c(4,3,1,2),oma=c(0,3,0,0),mgp=c(2,1,0))
  for (i in 1:length(vars)) {
i <<- i
print((vars)[i])
partialPlot(sel.rf,xdata,vars[i],which.class=1,xlab=vars[i],main="")
mtext("(Logit of probability of high severity)/2", side=2, line=1,
outer=T)
  }
}

Thanks again for the help.
-Greg Dillon




"Tim Howard"  
11/16/2010 05:36 AM

To
, 
cc

Subject
Re: [R] Force evaluation of variable when calling partialPlot






Greg, 
Two thoughts: 
1. It might be possible that 'vars' is a reserved word of sorts and if you 
change the name of your vector RF might be happier
2. A way that works for me is to call importance as follows:
 
sel.imp <- importance(sel.rf, class=NULL, scale=TRUE, type=NULL)
 
and then use the 'names' of the imp data frame to be absolutely clear to 
RF you are talking about the same variables
 
for(i in length(sel.imp){
 
partialPlot(sel.rf,xdata,names(sel.imp[i]),which.class=1,xlab=vars[i],main="")
}
 
 
Hope that helps. 
Tim Howard
 
>>>
Date: Mon, 15 Nov 2010 12:29:08 -0800 (PST)
From: gdillon 

RE: the folloing original question:
> I'm using the randomForest package and would like to generate partial
> dependence plots, one after another, for a variety of variables:
> 
> m <- randomForest( s, ... )
> varnames <- c( "var1", "var2", "var3", "var4" )   # var1..4 are all in
> data frame s
> for( v in varnames ) {
>partialPlot( x=m, pred.data=s, x.var=v )
> }
> 
> ...but this doesn't work, with partialPlot complaining that it can't
> find the variable "v".

I'm having a very similar problem using the partialPlot function. A
simplified version of my code looks like this:

data.in <- paste(basedir,"/sw_climate_dataframe_",root.name,".csv",sep="")
data <- read.csv(data.in)
vars <- c("sm1","precip.spring","tmax.fall","precip.fall")  #selected
variables in data frame
xdata <- as.data.frame(data[,vars])
ydata <- data[,5]
ntree <- 2000

rf.pdplots <- function() {
  sel.rf <- randomForest(xdata,ydata,ntree=ntree,keep.forest=TRUE)
 
par(family="sans",mfrow=c(2,2),mar=c(4,3,1,2),oma=c(0,3,0,0),mgp=c(2,1,0))
  for (i in 1:length(vars)) {
print((vars)[i])
partialPlot(sel.rf,xdata,vars[i],which.class=1,xlab=vars[i],main="")
mtext("(Logit of probability of high severity)/2", side=2, line=1,
outer=T)
  }
}

rf.pdplots()

When I run this code, with partialPlots embedded in a function, I get the
following error:
Error in eval(expr, envir, enclos) : object "i" not found

If I just take the code inside the function and run it (not embedded in 
the
function), it runs just fine. Is there some reason why partialPlots 
doesn't
like to be called from inside a function?

Other things I've tried/noticed:
1. If I comment out the line with the partialPlots call (and the next line
with mtext), the function runs as expected and prints the variable names 
one
at a time.
2. If the variable i is defined as a number (e.g., 4) in the global
environment, then the function will run, the names print out one at a 
time,
and four plots are created. HOWEVER, the plots are all for the last (4th)
variable, BUT the x labels actually are different on each plot (i.e., the
xlab is actually looping through the four values in vars).

Can anyone help me make sense of this? Thanks.

-Greg
-- 

[[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 vs print ??

[Joshua Wiley]

Hi,

It really depends, but almost always to plot a graph you would use
plot() (or some similar graphing function).  print() is usually more
for displaying data in the console, not graphically.  There are some
cases (e.g., inside a for loop), where you need to wrap the call to
plot (or similar graphing function) in print() (lattice graphics is
probably the most common example here).  Here is an
example/demonstration of plot() and print() in the context of a
histogram (which was calculated and saved into an object called 'x'):

## do histogram calculations but do not plot
x <- hist(rnorm(1000), plot = FALSE)

## print() just shows you the data
print(x)
## put data in a graph
plot(x)

## but you could also just use hist() without plot() or print()
hist(runif(100))

HTH,

Josh

On Tue, Nov 16, 2010 at 9:09 AM, skan  wrote:
>
> Hello
>
> What's the differente betwen using "plot" and using "print"  in order to
> plot a graph?
> For example in order to plot the result of a histogram.
>
> cheers
> --
> View this message in context: 
> http://r.789695.n4.nabble.com/plot-vs-print-tp3045256p3045256.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.
>



-- 
Joshua Wiley
Ph.D. Student, Health Psychology
University of California, Los Angeles
http://www.joshuawiley.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] R article in forbes/ somewhat OT

[Gabor Grothendieck]

On Tue, Nov 16, 2010 at 11:20 AM, Mike Marchywka  wrote:
> In many cases, the problem is not the R code but actually
> getting the data you need to make real-life decisions or form opinions.
> I sent a link to one of my rants to the OP on this topic
> but it is a bit political so I will spare the whole list.
> I would mention however that many day-to-day uses of R
> can be limited due to lack of good data in a computer
> readable format- you see people posting here about scraping
> stock quotes from human readable web pages but sometimes you need to do that 
> with
> data from government sources or they just cater to commercial
> products. If you use R at home, you may want to take note
> of requests from data sources on ideas about how they can optimize their
> data delivery for your needs.

You might want to check out the F. Scott Fitzgerald quote at the sqldf
home page.



-- 
Statistics & Software Consulting
GKX Group, GKX Associates Inc.
tel: 1-877-GKX-GROUP
email: ggrothendieck at gmail.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] plot vs print ??

[Douglas Bates]

On Tue, Nov 16, 2010 at 11:09 AM, skan  wrote:

> What's the differente betwen using "plot" and using "print"  in order to
> plot a graph?
> For example in order to plot the result of a histogram.

Could you give an example?

It seems that you are referring to graphics functions in packages such
as lattice and ggplot2 that produce objects that must be print'ed or
plot'ed before they are rendered on a graphics device.  However, we
would need to know if you mean the result of histogram() in the
lattice package or hist() in the graphics package before we could
answer.

__
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 Read a Large CSV into a Database with R

[Seth Falcon]

Hi Abhijit,

[I've cc'd R-help to keep the discussion on the list]

On Tue, Nov 16, 2010 at 8:06 AM, Abhijit Dasgupta
 wrote:
> Seth,
>
> I was looking for something like this too. I've a question. If
> you're reading the data from a connection, does R start reading the
> next chunk of data right after the previous chunk, or do we need to
> keep track of things using "skip"

The purpose of using a file connection is to allow R to keep its place
in the file as it reads and not have to re-read or skip.  This is
considerably more efficient.



-- 
Seth Falcon | @sfalcon | http://userprimary.net/

__
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 vs print ??

[David Winsemius]

http://cran.r-project.org/doc/FAQ/R-FAQ.html#Why-do-lattice_002ftrellis-graphics-not-work_003f

On Nov 16, 2010, at 12:09 PM, skan wrote:



Hello

What's the differente betwen using "plot" and using "print"  in  
order to

plot a graph?
For example in order to plot the result of a histogram.

cheers
--



David Winsemius, MD
West Hartford, CT

__
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] Non-positive definite cross-covariance matrices

[Jeff Bassett]

Giovanni,

Both matrices describing the points (A and B in my example) are the
same size, so the resulting matrix will always be square.  Also, the
equation I'm using is essentially the following identity:

Var(A + B) = Var(A) + Var(B) + Cov(A, B) + Cov(B, A)

All the covariance matrices that result from the Var() terms should be
positive definite, and while it seems possible that either of those
resulting from the Cov() terms may not be, the sum of the two should.
Do you agree?

Of course, the above identity only holds if the data is normally
distributed.  The Mardia test for multivariate normality in fact shows
that my data is not.  This may ultimately be my problem.  So maybe I
should be asking if you could point me toward some packages that can
transform my data so that it is normally distributed.

- Jeff


On Mon, Nov 15, 2010 at 3:56 PM, Giovanni Petris  wrote:
> What made you think that a cross-covariance matrix should be positive
> definite? Id does not even need to be a square matrix, or symmetric.
>
> Giovanni Petris
>
> On Mon, 2010-11-15 at 12:58 -0500, Jeff Bassett wrote:
>> I am creating covariance matrices from sets of points, and I am having
>> frequent problems where I create matrices that are non-positive
>> definite.  I've started using the corpcor package, which was
>> specifically designed to address these types of problems.  It has
>> solved many of my problems, but I still have one left.
>>
>> One of the matrices I need to calculate is a cross-covariance matrix.
>> In other words, I need to calculate cov(A, B), where A and B are each
>> a matrix defining a set of points.  The corpcor package does not seem
>> to be able to perform this operation.
>>
>> Can anyone suggest a way to create cross-covariance matrices that are
>> guaranteed (or at least likely) to be positive definite, either using
>> corpcor or another package?
>>
>> I'm using R 2.8.1 and corpcor 1.5.2 on Mac OS X 10.5.8.
>>
>> - Jeff
>>
>> __
>> 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] glmer, Error: Downdated X'X is not positive definite,49

[Douglas Bates]

It is often more effective to send questions about lmer or glmer to
the r-sig-mixed-mod...@r-project.org mailing list, which I am cc:ing
on this response.

On Tue, Nov 16, 2010 at 3:25 AM, Annika  wrote:
>
> Dear list,
>
> I am new to this list and I am new to the world of R. Additionally I am not
> very firm in statistics either but have to deal. So here is my problem:
> I have a dataset (which I attach at the end of the post) with a binomial
> response variable (alive or not) and three fixed factors
> (trapping,treat,sex). I do have repeated measures and would like to include
> one (enclosure) random factor.
> I chose to use the glmer function (and tried the lmer as well).
>
> m1<-glmer(alive~treat*sex*trapping*ID+(1|enclosure),family=binomial,data=survadult)

Is ID a factor?  In the data you posted it is an integer value but I
am guessing it should be a factor.  If so, you can't incorporate it in
the model because there is only one binary observation per level of
ID.

If I omit the ID in the fixed-effects part of the formula I can fit
the model with glmer.

> summary(m1 <- glmer(alive ~ trapping * treat * sex + (1|enclosure), trap, 
> binomial))
Generalized linear mixed model fit by maximum likelihood ['summary.mer']
 Family: binomial
Formula: alive ~ trapping * treat * sex + (1 | enclosure)
   Data: trap
 AIC  BIC   logLik deviance
213.2111 334.5440 -73.6056 147.2111

Random effects:
 GroupsNameVariance Std.Dev.
 enclosure (Intercept) 1.3721.171
Number of obs: 292, groups: enclosure, 13

Fixed effects:
 Estimate Std. Error z value
(Intercept)  10.61897   45.92452   0.231
trappingtrapping2 0.03992   65.59766   0.001
trappingtrapping3-0.01196   64.74665   0.000
trappingtrapping4-0.74441   71.77010  -0.010
treatm/p -0.26605   72.78328  -0.004
treatp/m -0.18007   81.63861  -0.002
treatp/p  0.50290   64.45725   0.008
sexm -8.61299   45.92746  -0.188
trappingtrapping2:treatm/p   -0.03657  103.39000   0.000
trappingtrapping3:treatm/p   -0.04916  101.87680   0.000
trappingtrapping4:treatm/p  -10.26432   91.32210  -0.112
trappingtrapping2:treatp/m   -0.04047  115.80425   0.000
trappingtrapping3:treatp/m   -8.59631   93.53334  -0.092
trappingtrapping4:treatp/m   -7.92162   98.52543  -0.080
trappingtrapping2:treatp/p   -0.05942   91.39839  -0.001
trappingtrapping3:treatp/p   -7.83591   78.98418  -0.099
trappingtrapping4:treatp/p   -9.60581   84.83662  -0.113
trappingtrapping2:sexm0.92624   65.61333   0.014
trappingtrapping3:sexm   -1.45692   64.75568  -0.022
trappingtrapping4:sexm   -0.88035   71.77685  -0.012
treatm/p:sexm 8.69446   93.57164   0.093
treatp/m:sexm 8.61957  106.06371   0.081
treatp/p:sexm 0.86615   64.46538   0.013
trappingtrapping2:treatm/p:sexm  -9.03327  118.95919  -0.076
trappingtrapping3:treatm/p:sexm  -7.51474  117.64003  -0.064
trappingtrapping4:treatm/p:sexm   0.81873  108.62671   0.008
trappingtrapping2:treatp/m:sexm  -9.51400  134.16089  -0.071
trappingtrapping3:treatp/m:sexm   0.48425  115.47998   0.004
trappingtrapping4:treatp/m:sexm  -2.37780  119.56181  -0.020
trappingtrapping2:treatp/p:sexm  -0.93097   91.42477  -0.010
trappingtrapping3:treatp/p:sexm   7.55060   79.00487   0.096
trappingtrapping4:treatp/p:sexm   9.47461   84.85476   0.112

> Both give the following error message:
>       Error in mer_finalize(ans) : Downdated X'X is not positive definite,
> 49.
>       Additionally: Warning:
>       glm.fit: fitted probabilities numerically 0 or 1 occurred
> I get that there is something wrong in my dataframe, but I don`t understand
> what and how I could change it.
> If I run a glm without the random factor
>       m1<-glm(alive~treat*sex*trapping*ID,family=binomial,data=survadult)
> it works out really fine, but is the wrong analysis in my opinion. (Does glm
> actually consider the repeated measures? Does glmer?)
>
> Is glmer the right function for me and how do I deal with the error I get?
>
> Thanks for any help for a non-professional :)
>
> Annika
>
>
> trapping        treat   sex     alive   ID      enclosure
>        trapping1       m/m     f       1       400     1
>        trapping1       m/m     f       1       401     1
>        trapping1       m/m     f       1       402     1
>        trapping1       m/m     m       0       403     1
>        trapping1       m/m     m       1       404     1
>        trapping1       m/m     m       1       405     1
>        trapping1       m/p     f       1       406     2
>        trapping1       m/p     f       1       407     2
>        trapping1       m/p     f       1       408     2
>        trapping1       m/p     m       1       409     2
>        trapping1       m/p     m       1      

Re: [R] Simple error handling in R

[Gabor Grothendieck]

On Tue, Nov 16, 2010 at 10:58 AM, Aleksey Naumov  wrote:
> Hi R experts,
>
> I am looking for a simple error handling approach, whereby I could stop
> function execution with a customized error message. For example:
>
> for (i in 1:10) {
>   if (i == 5)
>      # I'd like to be able to stop right here with an error message I have
> complete control over
> }
>
> The problem with stop() is that I cannot control completely what gets
> printed to the terminal, even with stop(call.=FALSE) there is the "Error:"
> string. I've worked through examples (whatever few there are) for try() and
> tryCatch() and I still cannot understand how to do this. If I supply my own
> error handler function with tryCatch(..., error=function(e) ...) I can
> control the error message, but the loop continues on to i=6, etc.
>
> So I am struggling with error handling in R... It' seems its a lot simpler
> and more consistent e.g. in Python.
>
> Any help would be greatly appreciated!
>

There was a long thread on this a couple of years ago.   You might
want to look at this response regarding callCC and also the other
responses:

http://tolstoy.newcastle.edu.au/R/e5/help/08/11/6967.html

-- 
Statistics & Software Consulting
GKX Group, GKX Associates Inc.
tel: 1-877-GKX-GROUP
email: ggrothendieck at gmail.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] rotate column names in large matrix

[Chris Stubben]

Just increase the margins on the left side and add the rownames

x <- cor(matrix(rnorm(600), 60, 100))
rownames(x)<-paste("row", 1:100)
op<-par(mar=c(1,5,1,1), xpd=TRUE)
image(t(x[nrow(x):1,]), axes=FALSE)
 text(-0.01, seq(0,1,length=nrow(x) ), rownames(x), pos = 2,  offset = 0,
cex = .7)

Another option is to search the R graphics manual at
http://rgm2.lab.nig.ac.jp/RGM2/images.php and maybe you'll find a function
in all those packages that plots what you need.   Try searching for "symnum"
or "image key" or "matrix plot" (in quotes).  For "matrix plot", there's a
function in the plotrix package that will add a key or even show values in
the "cells"..

library(plotrix)
color2D.matplot(x, show.legend=TRUE, red=1, blue=0)

Chris




Lara Poplarski wrote:
> 
> Chris, could you please suggest how to modify what you sent to also show
> the
> same labels as (horizontal) row names? I have not yet mastered the details
> of R graphics...
> 

-- 
View this message in context: 
http://r.789695.n4.nabble.com/rotate-column-names-in-large-matrix-tp3043493p3045366.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] help using function merge_all()

[Joshua Wiley]

Hi,

Did you also load the package (you'll need to load it every session
you want to use it)?

## load reshape package
library(reshape)
## now try
merge_all(Br, Ki, Lu, Pr, by="Genes")


If this does not resolve your problem, please run:   sessionInfo()
at your console and report the output (it will tell us information
about your version of R, loaded packages, etc.).

Cheers,

Josh

On Tue, Nov 16, 2010 at 8:33 AM, arjun  wrote:
>
> Hi,
>       I want to merge 4 data frames with one column in common but I am
> getting error message while using this function. Can any one help me out.
>> merge_all(Br,Ki,Lu,Pr,by="Genes")
> Error: could not find function "merge_all"
>
> I have installed the package: reshape but I still get this error
> --
> View this message in context: 
> http://r.789695.n4.nabble.com/help-using-function-merge-all-tp3045173p3045173.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.
>



-- 
Joshua Wiley
Ph.D. Student, Health Psychology
University of California, Los Angeles
http://www.joshuawiley.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] Please help me on a simple avg question

[facehappywy]

Thanks every one first.

I have a dataset like this

YearMonth   Day HourDIR SPD (m/s)   SKYCOVERTEMP 
(C)RH(%)
20091   1   0   310 13.858243   -5  23
20091   1   1   330 14.305283   -7.2
26
20091   1   2   320 10.728963   -7.2
24
20091   1   3   320 10.728963   -6.7
23
20091   1   4   320 10.281923   -7.2
26
20091   1   5   330 10.281923   -7.2
26
20091   1   6   330 10.281923   -7.2
25
20091   1   6   320 12.964163   -7.2
26
20091   1   7   330 11.176  3   -7.226
20091   1   8   330 10.728963   -7.8
28
20091   1   9   330 10.728963   -7.8
27
20091   1   9   320 11.176  3   -7.226
20091   1   10  320 12.517123   -7.8
24
20091   1   11  310 9.83488 3   -7.824
20091   1   12  310 9.83488 3   -8.325
20091   1   12  310 10.728963   -8.9
26
20091   1   13  320 9.83488 3   -7.826
20091   1   14  310 9.83488 3   -7.223
20091   1   15  310 9.83488 3   -7.223
20091   1   15  320 10.728960   -7.2
23
20091   1   16  320 8.04672 0   -6.121
20091   1   17  300 9.83488 0   -3.918
20091   1   18  300 9.83488 0   -4.419
20091   1   18  310 9.38784 0   -2.816
20091   1   19  320 7.59968 3   -2.816
20091   1   20  310 6.25856 3   -2.818
20091   1   21  310 6.25856 3   -2.818
20091   1   21  310 7.15264 3   -2.816
20091   1   22  310 7.59968 3   -2.816
20091   1   23  300 5.81152 3   -3.916


It is a one year (2009), 12 month, every day, every hour data. Notice that
there are sometimes several readings  within one hour. I want to calculate
the hourly avg for each day in each month. Anyone could help me on this. 
-- 
View this message in context: 
http://r.789695.n4.nabble.com/Please-help-me-on-a-simple-avg-question-tp3045440p3045440.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] Non-positive definite cross-covariance matrices

[Peter Langfelder]

On Tue, Nov 16, 2010 at 9:40 AM, Jeff Bassett  wrote:
> Giovanni,
>
> Both matrices describing the points (A and B in my example) are the
> same size, so the resulting matrix will always be square.  Also, the
> equation I'm using is essentially the following identity:
>
> Var(A + B) = Var(A) + Var(B) + Cov(A, B) + Cov(B, A)
>
> All the covariance matrices that result from the Var() terms should be
> positive definite, and while it seems possible that either of those
> resulting from the Cov() terms may not be, the sum of the two should.
> Do you agree?
>
> Of course, the above identity only holds if the data is normally
> distributed.

Hi Jeff,

as far as I can see, the identity above is an identity and holds
irrespective of how your data is distributed (just write out the
difinition of var and cov and you will see that the equation is always
true).

It is easy to come up with examples where Cov(A, B) + Cov(B, A) is not
positive definite. As an extreme example, consider a matrix A (say 10
columns, 100 rows) such that the off-diagonal covariances are all zero
and the columns are standardized, so cov(A) = diag(1, 1, 1, ...). Then
take B = -A, so cov(A, B) = cov(B, A) = diag(-1, -1, -1, ...).
Obviously, cov(A, B) + cov(B, A) is not positively definite, in fact
it is negative definite.

If the matrices A and B are completely independent, it is not very
likely that cov(A,B) + cov(B,A) will be positive definite. For
example, if A and B had just one column, there's a 50-50 chance that
cov(A, B) is negative (single number). When you have more than one
column, the chance is even less than 50-50 because all eigenvalues
would have to be positive.

You may be able to generate matrices whose cov(A, B) is positive
definite by starting with a matrix A, then generating a random matrix
B0, subtracting from the columns of B0 the projections of columns of
B0 into the columns of A, then adding a small multiple of A to get B.
But this may not be what you want or need.

Alternatively (and more likely), something in your approach may need
re-thinking.

HTH,

Peter


> - Jeff
>
>
> On Mon, Nov 15, 2010 at 3:56 PM, Giovanni Petris  wrote:
>> What made you think that a cross-covariance matrix should be positive
>> definite? Id does not even need to be a square matrix, or symmetric.
>>
>> Giovanni Petris
>>
>> On Mon, 2010-11-15 at 12:58 -0500, Jeff Bassett wrote:
>>> I am creating covariance matrices from sets of points, and I am having
>>> frequent problems where I create matrices that are non-positive
>>> definite.  I've started using the corpcor package, which was
>>> specifically designed to address these types of problems.  It has
>>> solved many of my problems, but I still have one left.
>>>
>>> One of the matrices I need to calculate is a cross-covariance matrix.
>>> In other words, I need to calculate cov(A, B), where A and B are each
>>> a matrix defining a set of points.  The corpcor package does not seem
>>> to be able to perform this operation.
>>>
>>> Can anyone suggest a way to create cross-covariance matrices that are
>>> guaranteed (or at least likely) to be positive definite, either using
>>> corpcor or another package?
>>>
>>> I'm using R 2.8.1 and corpcor 1.5.2 on Mac OS X 10.5.8.
>>>
>>> - Jeff
>>>
>>> __
>>> 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] Please help me on a simple avg question

[Henrique Dallazuanna]

Try this:

aggregate(TEMP ~ Year + Month + Day, DF, mean)

On Tue, Nov 16, 2010 at 4:48 PM, facehappywy  wrote:

>
> Thanks every one first.
>
> I have a dataset like this
>
> YearMonth   Day HourDIR SPD (m/s)   SKYCOVER
>  TEMP (C)RH(%)
> 20091   1   0   310 13.858243   -5  23
> 20091   1   1   330 14.305283
> -7.226
> 20091   1   2   320 10.728963
> -7.224
> 20091   1   3   320 10.728963
> -6.723
> 20091   1   4   320 10.281923
> -7.226
> 20091   1   5   330 10.281923
> -7.226
> 20091   1   6   330 10.281923
> -7.225
> 20091   1   6   320 12.964163
> -7.226
> 20091   1   7   330 11.176  3   -7.226
> 20091   1   8   330 10.728963
> -7.828
> 20091   1   9   330 10.728963
> -7.827
> 20091   1   9   320 11.176  3   -7.226
> 20091   1   10  320 12.517123
> -7.824
> 20091   1   11  310 9.83488 3   -7.824
> 20091   1   12  310 9.83488 3   -8.325
> 20091   1   12  310 10.728963
> -8.926
> 20091   1   13  320 9.83488 3   -7.826
> 20091   1   14  310 9.83488 3   -7.223
> 20091   1   15  310 9.83488 3   -7.223
> 20091   1   15  320 10.728960
> -7.223
> 20091   1   16  320 8.04672 0   -6.121
> 20091   1   17  300 9.83488 0   -3.918
> 20091   1   18  300 9.83488 0   -4.419
> 20091   1   18  310 9.38784 0   -2.816
> 20091   1   19  320 7.59968 3   -2.816
> 20091   1   20  310 6.25856 3   -2.818
> 20091   1   21  310 6.25856 3   -2.818
> 20091   1   21  310 7.15264 3   -2.816
> 20091   1   22  310 7.59968 3   -2.816
> 20091   1   23  300 5.81152 3   -3.916
>
>
> It is a one year (2009), 12 month, every day, every hour data. Notice that
> there are sometimes several readings  within one hour. I want to calculate
> the hourly avg for each day in each month. Anyone could help me on this.
> --
> View this message in context:
> http://r.789695.n4.nabble.com/Please-help-me-on-a-simple-avg-question-tp3045440p3045440.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.
>



-- 
Henrique Dallazuanna
Curitiba-Paraná-Brasil
25° 25' 40" S 49° 16' 22" O

[[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] Error: invalid type(list) for variable when using lm()

[Tianchan Niu]

Dear All,I would like to do multiple regression in R. I used: lm(y~X),
where y is a n by 1 vector, and X is a n by m matrix. I kept getting the error
message:Error in model.frame.default(formula = y ~ X,  :  invalid type (list) 
for variable 'X'. However, when I used: lm(y~X[,1]+X[,2]+X[,3]+…+X[,m]), it
works well, but this is not the form I prefer, it makes my codes 
complicated.Please help.Thank you very much,Jane__
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] Pass character vector to function argument

[David Winsemius]

On Nov 16, 2010, at 1:04 PM, Kevin Ummel wrote:

Sorry, I shouldn't have used 'sum' as an example; I am looking for a  
solution in the case of functions that do not result in simple  
vectors or matrices.


The "real-world" example is something like the following using the  
'raster' package, though (I think) any command producing an S4  
object could be a candidate:


r1 <- raster(system.file("external/test.grd", package="raster"))*1
r2 <- raster(system.file("external/test.grd", package="raster"))*2

obs=objects(pattern=glob2rx("r?"))

#Would like to generate the equivalent of...

stack(r1,r2)


Have you tried what would seem to be the obvious generalization of the  
example offered:


str(stack(sapply(obs, get)))   # it does not throw any error and seems  
to have data in it




...using a variation of obs as an argument to 'stack' (potentially  
very long).


Also tested:

str(
stack( c(sapply(obs, get),sapply(obs, get)) )  # doing a four object  
test

)  # no error

(I did not see any indication in the help page for get that it's  
operation was restricted to S3 objects. If there is such a restriction  
then warnings should be added there.)





Cheers,
Kevin


On Nov 16, 2010, at 5:07 PM, David Winsemius wrote:



On Nov 16, 2010, at 10:34 AM, Kevin Ummel wrote:

A bit embarrassed to post this seemingly trivial question, but I  
can't find anything in the archive that's quite relevant:


a1=1
a2=2

obs=objects(pattern=glob2rx("a?"))

I want to utilize 'obs' as a function argument to produce  
something like:


sum(a1,a2)



?get

sum(sapply(obs, get))

[1] 3


Obviously, sum(obs) doesn't work, but I've tried variations of  
'eval', 'parse', 'mget', and 'noquote' without success. What am I  
missing?


Thanks,
Kevin


--
David Winsemius, MD
West Hartford, CT





David Winsemius, MD
West Hartford, CT

__
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] DBLEPR?

[Prof. John C Nash]

Ravi Varadhan and I have been looking at UCMINF to try to identify why it gives 
occasional
(but not reproducible) errors, seemingly on Windows only. There is some 
suspicion that its
use of DBLEPR for finessing the Fortran WRITE() statements may be to blame. 
While I can
find DBLEPR in Venables and Ripley, it doesn't get much mention after about 
2000 in the
archives, though it is in the R FAQ and Brian R. mentions they are in libR in a 
2009 post.
Are the xxxPR routines now deprecated (particularly for 64 bit systems) or 
still OK to
use?  If OK, can anyone point to documentation and examples?

JN

__
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: invalid type(list) for variable when using lm()

[Henrique Dallazuanna]

Try this:

y <- rnorm(100)
X <- matrix(runif(100 * 10), ncol = 10)

lm(y ~ ., data = cbind.data.frame(y, X))


On Tue, Nov 16, 2010 at 5:07 PM, Tianchan Niu wrote:

> Dear All,I would like to do multiple regression in R. I used: lm(y~X),
> where y is a n by 1 vector, and X is a n by m matrix. I kept getting the
> error
> message:Error in model.frame.default(formula = y ~ X,  :  invalid type
> (list) for variable 'X'. However, when I used:
> lm(y~X[,1]+X[,2]+X[,3]+…+X[,m]), it
> works well, but this is not the form I prefer, it makes my codes
> complicated.Please help.Thank you very much,Jane
> __
> 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.
>
>


-- 
Henrique Dallazuanna
Curitiba-Paraná-Brasil
25° 25' 40" S 49° 16' 22" O

[[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] Error: invalid type(list) for variable when using lm()

[Phil Spector]

Tianchan -
   Your X is not a matrix -- it's a dataframe.  Probably the 
simplest solution is to use


lm(y~as.matrix(X))

but you should also learn the difference between a data frame
and a matrix.

- Phil Spector
 Statistical Computing Facility
 Department of Statistics
 UC Berkeley
 spec...@stat.berkeley.edu


On Tue, 16 Nov 2010, Tianchan Niu wrote:


Dear All,I would like to do multiple regression in R. I used: lm(y~X),
where y is a n by 1 vector, and X is a n by m matrix. I kept getting the error
message:Error in model.frame.default(formula = y ~ X,  :  invalid type (list) 
for variable 'X'. However, when I used: lm(y~X[,1]+X[,2]+X[,3]+…+X[,m]), it
works well, but this is not the form I prefer, it makes my codes 
complicated.Please help.Thank you very much,Jane__
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] Extending a plot in a loop (With attachment)

[Nasrin Pak]

My problem is that I have a data set for every day of measurement in a
seperate file and I want to plot one parameter of the data for all the days
in one graph. I tried to use for loop but only the last data remains in the
program memory, I don`t know how to plot each day`s data continusly after
the others(or how to extending the x axis.) Would you please help me with
it?

This a plot for one day:

radiation.data
<-read.table("C:/updated_CFL_Rad_files/2008/RAD_2008_JD101_0410.dat", header
= TRUE,sep = ",", quote = " ", dec = ".")
> attach(radiation.data)
The following object(s) are masked from 'radiation.data (position 3)':

Batt_avg, Batt_st, Day, Hour, Kdown_avg, Kdown_st, LW_in, LW_in_st,
Minute, Month, PanelT_avg, PanelT_st, PAR_avg, PAR_st, Sec,
Tcase_avg, Tcase_st, Tdome_avg, Tdome_st, Thermopile_avg,
Thermopile_st, Tuv_avg, Tuv_st, Uva_avg, Uva_st, Uvb_avg, Uvb_st,
Year
> names(radiation.data)
 [1] "Year"   "Month"  "Day""Hour"
 [5] "Minute" "Sec""Batt_avg"   "PanelT_avg"
 [9] "Batt_st""PanelT_st"  "Kdown_avg"  "Thermopile_avg"
[13] "Tcase_avg"  "Tdome_avg"  "LW_in"  "PAR_avg"
[17] "Tuv_avg""Uvb_avg""Uva_avg""Kdown_st"
[21] "Thermopile_st"  "Tcase_st"   "Tdome_st"   "LW_in_st"
[25] "PAR_st" "Tuv_st" "Uvb_st" "Uva_st"

 plot(((PAR_avg*0.216)/Uvb_avg),
main="Par/UVB",xlab="minutes",ylab="Par/UVB")


and this is the algorithm I tried  for plotting all the data in one plot:

x<- matrix( list.files("C:/updated_CFL_Rad_files",full=TRUE)) # putting all
data sets in a matrix
  for(i in 1:100) {
  if(i < 101) next
 radiation.data <-read.table(x[i], header = TRUE,sep = ",", quote = " ",
dec = ".")
 attach(radiation.data)
 plot(i*Hour*60+Minute,PAR_avg,main="PAR",xlab="Hour",ylab="Par")
dev.print(device=postscript, "C:/graph5.eps", onefile=FALSE,
horizontal=FALSE)
   }
The plot I see is the last file's plot, I don't know how to keep previous
data and continue within the same plot.

* I have attached a sample of the data to this email.

-- 
Sincerely

Nasrin  Pak







-- 
Sincerely

Nasrin  Pak
__
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] Pass character vector to function argument

[Kevin Ummel]

Sorry, I shouldn't have used 'sum' as an example; I am looking for a solution 
in the case of functions that do not result in simple vectors or matrices.

The "real-world" example is something like the following using the 'raster' 
package, though (I think) any command producing an S4 object could be a 
candidate:

r1 <- raster(system.file("external/test.grd", package="raster"))*1
r2 <- raster(system.file("external/test.grd", package="raster"))*2

obs=objects(pattern=glob2rx("r?"))

#Would like to generate the equivalent of...

stack(r1,r2)

...using a variation of obs as an argument to 'stack' (potentially very long).

Cheers,
Kevin 


On Nov 16, 2010, at 5:07 PM, David Winsemius wrote:

> 
> On Nov 16, 2010, at 10:34 AM, Kevin Ummel wrote:
> 
>> A bit embarrassed to post this seemingly trivial question, but I can't find 
>> anything in the archive that's quite relevant:
>> 
>> a1=1
>> a2=2
>> 
>> obs=objects(pattern=glob2rx("a?"))
>> 
>> I want to utilize 'obs' as a function argument to produce something like:
>> 
>> sum(a1,a2)
>> 
> 
> ?get
> > sum(sapply(obs, get))
> [1] 3
> 
> 
>> Obviously, sum(obs) doesn't work, but I've tried variations of 'eval', 
>> 'parse', 'mget', and 'noquote' without success. What am I missing?
>> 
>> Thanks,
>> Kevin
>> 
> -- 
> David Winsemius, MD
> West Hartford, CT
> 

__
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] Anyone can help with this question

[Mariana]

Hi there:
I am a total beginner in R, and I have a simple question:
I have a table with thousands of lines that represent locations, and two
columns: latitude and longitude. I need to randomly sample 1000 lines. How
do I do it? I know the command "sample", but it samples elements
independently, not lines. 
If there is a better place for me to ask that type of question, please let
me know.
Thanks a lot,
Mariana (from Brazil)

-- 
View this message in context: 
http://r.789695.n4.nabble.com/Anyone-can-help-with-this-question-tp3045265p3045265.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] discriminant function analysis

[Mike Gibson]

My objective is to look at differences in two species of fish from morphometric 
measurements.  My morphometric measurements are head length, eye diameter, 
snout length, and measurements from tail to each fin. I want to use 
discrimanant function analyis to determine if there are differences between the 
two species.  
 
I am familiar with R but new to discrimannt function analysis.  I want to learn 
how to run this analysis in R.  My internet search has not been very 
successful.  I did see refrence of the linear discriminant function (lda) in R 
but I am not sure if this is the correct function for me.  
 
Any advice and especially a reference to any examples would be greatly 
appreciated.  
 
Thanks.  
 
Mike  
[[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] Odp: Sampling problem

[wangwallace]

Fabulicious! It worked!!! 

One more question, in the following data frame as posted above:

SubIDCSE1 CSE2 CSE3 CSE4 WSE1 WSE2 WSE3 WSE4
  1  6  5   6   2  6  22   4
  2  6  4   7   2  6  62   3
  3  5  5   5   5  5  54   5
  4  5  4   3   4  4  45   2
  5  5  6   7   5  6  44   1
  6  5  4   3   6  4  37   3
  7  3  6   6   3  6  52   1
  8  3  6   6   3  6  54   7 

I want to draw the first random sample consisting of a row of integers under
the first group of variables (CSE1, CSE2, CSE3, CSE4). For example, assuming
the first draw yielded a sample of the first row (6, 5, 6, 2), now I want to
draw another random sample consisting of two rows of integers under the
second group of variables (WSE1, WSE2, WSE3, WSE4). Also, for the second
draw, I want to restrict a vector I am going to sample from to only those
rows that are not correspond to SubID I have sampled. That is, I want to
sample two rows of integers under the second group of variables (WSE1, WSE2,
WSE3, WSE4) from rows 2, 3, 4, 5, 6, 7, and 8.

Also, I want to repeat this whole process (drawing 1 random row of integers
under the first group of variables first, AND then another two random rows
under the second group of variables) for 1000 times. Any ideas? would that
be possible to do it by just revising the syntax you wrote above? many
thanks!!!
-- 
View this message in context: 
http://r.789695.n4.nabble.com/Sampling-problem-tp3043804p3045352.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] Playwith-problem with loading

[Fencl, Martin]

Helllo,
I am having trouble with running the library Playwith in the R-2.12.0. running 
under 32bit Windows XP. After calling the library the error message "The 
procedure entry point gdk_cairo_reset_clip could not be located in the dynamic 
library libgdk-win32-2.0-0.dll." occurs and the R asks for instaling the GTK+ 
package. However the instalation of GTK+ thorugh the R fails. I have installed 
GTK+ by installing GIMP 2.0 version 2.6.10 and set the path of environment 
variables. But it still does not work. It seems to be that the problem is in 
gtk2 (which I also installed).
Thank you for any hints!


> library(playwith)
Loading required package: lattice
Loading required package: cairoDevice
Loading required package: gWidgetsRGtk2
Loading required package: gWidgets
Error in inDL(x, as.logical(local), as.logical(now), ...) : 
  unable to load shared object 
'C:/PROGRA~1/R/R-212~1.0/library/RGtk2/libs/i386/RGtk2.dll':
  LoadLibrary failure:  The specified procedure could not be found.


Failed to load RGtk2 dynamic library, attempting to install it.
trying URL 
'http://downloads.sourceforge.net/gtk-win/gtk2-runtime-2.22.0-2010-10-21-ash.exe?download'
Content type 'application/x-msdos-program' length 7820679 bytes (7.5 Mb)
opened URL
downloaded 7.5 Mb

Learn more about GTK+ at http://www.gtk.org
If the package still does not load, please ensure that GTK+ is installed and 
that it is on your PATH environment variable
IN ANY CASE, RESTART R BEFORE TRYING TO LOAD THE PACKAGE AGAIN
Error : .onAttach failed in attachNamespace() for 'gWidgetsRGtk2', details:
  call: .Call(name, ..., PACKAGE = PACKAGE)
  error: C symbol name "S_gtk_icon_factory_new" not in DLL for package "RGtk2"
Error: package 'gWidgetsRGtk2' could not be loaded
> 

__
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] Anyone can help with this question

[David Winsemius]

On Nov 16, 2010, at 12:13 PM, Mariana wrote:



Hi there:
I am a total beginner in R,


Time to read the Posting Guide. (Especially since Nabble appears to  
deficient in providing an introductory document, at least judging from  
the email behavior of some of its users.)



and I have a simple question:
I have a table with thousands of lines that represent locations, and  
two
columns: latitude and longitude. I need to randomly sample 1000  
lines. How

do I do it? I know the command "sample", but it samples elements
independently, not lines.


So use that vector as an index:

latlong.tbl[ sample(NROW(latlong.tbl), 1000) , ]

# should work for matrices or dataframes or table objects

If there is a better place for me to ask that type of question,  
please let

me know.
Thanks a lot,
Mariana (from Brazil)



--
David Winsemius, MD
West Hartford, CT

__
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] Pass character vector to function argument

[Kevin Ummel]

Thanks, David. That does, indeed, work. It didn't occur to me that a list would 
do the job as an argument.

Thanks for the fix!

kevin


On Nov 16, 2010, at 6:58 PM, David Winsemius wrote:

> 
> On Nov 16, 2010, at 1:04 PM, Kevin Ummel wrote:
> 
>> Sorry, I shouldn't have used 'sum' as an example; I am looking for a 
>> solution in the case of functions that do not result in simple vectors or 
>> matrices.
>> 
>> The "real-world" example is something like the following using the 'raster' 
>> package, though (I think) any command producing an S4 object could be a 
>> candidate:
>> 
>> r1 <- raster(system.file("external/test.grd", package="raster"))*1
>> r2 <- raster(system.file("external/test.grd", package="raster"))*2
>> 
>> obs=objects(pattern=glob2rx("r?"))
>> 
>> #Would like to generate the equivalent of...
>> 
>> stack(r1,r2)
> 
> Have you tried what would seem to be the obvious generalization of the 
> example offered:
> 
> str(stack(sapply(obs, get)))   # it does not throw any error and seems to 
> have data in it
> 
>> 
>> ...using a variation of obs as an argument to 'stack' (potentially very 
>> long).
> 
> Also tested:
> 
> str(
> stack( c(sapply(obs, get),sapply(obs, get)) )  # doing a four object test
> )  # no error
> 
> (I did not see any indication in the help page for get that it's operation 
> was restricted to S3 objects. If there is such a restriction then warnings 
> should be added there.)
> 
> 
>> 
>> Cheers,
>> Kevin
>> 
>> 
>> On Nov 16, 2010, at 5:07 PM, David Winsemius wrote:
>> 
>>> 
>>> On Nov 16, 2010, at 10:34 AM, Kevin Ummel wrote:
>>> 
 A bit embarrassed to post this seemingly trivial question, but I can't 
 find anything in the archive that's quite relevant:
 
 a1=1
 a2=2
 
 obs=objects(pattern=glob2rx("a?"))
 
 I want to utilize 'obs' as a function argument to produce something like:
 
 sum(a1,a2)
 
>>> 
>>> ?get
 sum(sapply(obs, get))
>>> [1] 3
>>> 
>>> 
 Obviously, sum(obs) doesn't work, but I've tried variations of 'eval', 
 'parse', 'mget', and 'noquote' without success. What am I missing?
 
 Thanks,
 Kevin
 
>>> -- 
>>> David Winsemius, MD
>>> West Hartford, CT
>>> 
>> 
> 
> David Winsemius, MD
> West Hartford, CT
> 

__
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] format secondary axis for dates

[Jannis]

Dear List,


this may be a Newbi question and may have been asked several times, but i am 
too stupid to find the posts.


I have a plot of values against POSIXct time steps. If I want to add a second x 
axis to the top margin of the plot, only numbers are at the tickmarks. Is there 
a straightforward way to specify the format to convert them to character 
representations (similar to the labels at the bottom)?

x<-as.POSIXct(1:1000*(60^2),origin='01-01-1970')
y=rnorm(1000)
plot(x,y)
axis(3)


Thanks for your help
Jannis

 



__
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] Pass character vector to function argument

[Henrique Dallazuanna]

Using David's examples:

do.call(stack, lapply(obs, get))


On Tue, Nov 16, 2010 at 4:04 PM, Kevin Ummel  wrote:

> Sorry, I shouldn't have used 'sum' as an example; I am looking for a
> solution in the case of functions that do not result in simple vectors or
> matrices.
>
> The "real-world" example is something like the following using the 'raster'
> package, though (I think) any command producing an S4 object could be a
> candidate:
>
> r1 <- raster(system.file("external/test.grd", package="raster"))*1
> r2 <- raster(system.file("external/test.grd", package="raster"))*2
>
> obs=objects(pattern=glob2rx("r?"))
>
> #Would like to generate the equivalent of...
>
> stack(r1,r2)
>
> ...using a variation of obs as an argument to 'stack' (potentially very
> long).
>
> Cheers,
> Kevin
>
>
> On Nov 16, 2010, at 5:07 PM, David Winsemius wrote:
>
> >
> > On Nov 16, 2010, at 10:34 AM, Kevin Ummel wrote:
> >
> >> A bit embarrassed to post this seemingly trivial question, but I can't
> find anything in the archive that's quite relevant:
> >>
> >> a1=1
> >> a2=2
> >>
> >> obs=objects(pattern=glob2rx("a?"))
> >>
> >> I want to utilize 'obs' as a function argument to produce something
> like:
> >>
> >> sum(a1,a2)
> >>
> >
> > ?get
> > > sum(sapply(obs, get))
> > [1] 3
> >
> >
> >> Obviously, sum(obs) doesn't work, but I've tried variations of 'eval',
> 'parse', 'mget', and 'noquote' without success. What am I missing?
> >>
> >> Thanks,
> >> Kevin
> >>
> > --
> > David Winsemius, MD
> > West Hartford, CT
> >
>
> __
> 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.
>



-- 
Henrique Dallazuanna
Curitiba-Paraná-Brasil
25° 25' 40" S 49° 16' 22" O

[[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] Pass character vector to function argument

[Kevin Ummel]

On second glance, while it works for 'stack', it doesn't appear to work for 
'overlay':

> overlay(sapply(obs,get))
Error in function (classes, fdef, mtable)  : 
  unable to find an inherited method for function "overlay", for signature 
"list", "missing"

At this point, this may be more of a 'raster' issue for the R-Sig-Geo list...

On Nov 16, 2010, at 7:09 PM, Kevin Ummel wrote:

> Thanks, David. That does, indeed, work. It didn't occur to me that a list 
> would do the job as an argument.
> 
> Thanks for the fix!
> 
> kevin
> 
> 
> On Nov 16, 2010, at 6:58 PM, David Winsemius wrote:
> 
>> 
>> On Nov 16, 2010, at 1:04 PM, Kevin Ummel wrote:
>> 
>>> Sorry, I shouldn't have used 'sum' as an example; I am looking for a 
>>> solution in the case of functions that do not result in simple vectors or 
>>> matrices.
>>> 
>>> The "real-world" example is something like the following using the 'raster' 
>>> package, though (I think) any command producing an S4 object could be a 
>>> candidate:
>>> 
>>> r1 <- raster(system.file("external/test.grd", package="raster"))*1
>>> r2 <- raster(system.file("external/test.grd", package="raster"))*2
>>> 
>>> obs=objects(pattern=glob2rx("r?"))
>>> 
>>> #Would like to generate the equivalent of...
>>> 
>>> stack(r1,r2)
>> 
>> Have you tried what would seem to be the obvious generalization of the 
>> example offered:
>> 
>> str(stack(sapply(obs, get)))   # it does not throw any error and seems to 
>> have data in it
>> 
>>> 
>>> ...using a variation of obs as an argument to 'stack' (potentially very 
>>> long).
>> 
>> Also tested:
>> 
>> str(
>> stack( c(sapply(obs, get),sapply(obs, get)) )  # doing a four object test
>> )  # no error
>> 
>> (I did not see any indication in the help page for get that it's operation 
>> was restricted to S3 objects. If there is such a restriction then warnings 
>> should be added there.)
>> 
>> 
>>> 
>>> Cheers,
>>> Kevin
>>> 
>>> 
>>> On Nov 16, 2010, at 5:07 PM, David Winsemius wrote:
>>> 
 
 On Nov 16, 2010, at 10:34 AM, Kevin Ummel wrote:
 
> A bit embarrassed to post this seemingly trivial question, but I can't 
> find anything in the archive that's quite relevant:
> 
> a1=1
> a2=2
> 
> obs=objects(pattern=glob2rx("a?"))
> 
> I want to utilize 'obs' as a function argument to produce something like:
> 
> sum(a1,a2)
> 
 
 ?get
> sum(sapply(obs, get))
 [1] 3
 
 
> Obviously, sum(obs) doesn't work, but I've tried variations of 'eval', 
> 'parse', 'mget', and 'noquote' without success. What am I missing?
> 
> Thanks,
> Kevin
> 
 -- 
 David Winsemius, MD
 West Hartford, CT
 
>>> 
>> 
>> David Winsemius, MD
>> West Hartford, CT
>> 
> 

__
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] Computing Rolling Average

[Paolo Rossi]

Hi,
Can anyone suggest a  clever way to compute a rolling weekly average of the
columns in a matrix?  The column bit is straightforward – use apply given a
function which does what you want on a column. With regard to a particular
column, the obvious way is to run a for loop indexing the last 7 days and
computing the average . I simply would like to know if there is a  better /
quicker way.



Code:
Given a,
> a= array(1:100, dim = c(17,2))
> a
  [,1] [,2]
 [1,]1   18
 [2,]2   19
 [3,]3   20
 [4,]4   21
 [5,]5   22
 [6,]6   23
 [7,]7   24
 [8,]8   25
 [9,]9   26
[10,]   10   27
[11,]   11   28
[12,]   12   29
[13,]   13   30
[14,]   14   31
[15,]   15   32
[16,]   16   33
[17,]   17   34
one needs to start computing the average from obs 7 s (at obs 7 you have a
full week to compute the average) and compute the rolling weekly average
from day 7 onwards
Results will look like b
  [,1] [,2]
 [1,]4   14
 [2,]5   21
 [3,]6   22
 [4,]7   23
 [5,]8   24
 [6,]9   25
 [7,]   10   26
 [8,]   11   27
 [9,]   12   28
[10,]   13   29
Thanks in advance,

Paolo



On 15 November 2010 20:59, wangwallace  wrote:

>
> Hey,
>
> I am hoping someone can help me with a sampling question.
>
> I have a data frame of 8 variables (the first column is the subjects' id):
>
>SubIDCSE1 CSE2 CSE3 CSE4 WSE1 WSE2 WSE3 WSE4
>  1  6  5   6   2  6  22   4
>  2  6  4   7   2  6  62   3
>  3  5  5   5   5  5  54   5
>  4  5  4   3   4  4  45   2
>  5  5  6   7   5  6  44   1
>  6  5  4   3   6  4  37   3
>  7  3  6   6   3  6  52   1
>  8  3  6   6   3  6  54   7
>
> the 6 variables are categorized into two groups with CSE1, CSE2, CSE3, and
> CSE4 in one group and the rest in another group.
>
> >sample(data[,2:4],2,replace=FALSE)
>
>   CSE1 CSE2
> 1  65
> 2  64
> 3  55
> 4  54
> 5  56
> 6  54
> 7  36
> 8  36
>
> Now I want to sample 1 column from another group of variables (i.e., WSE1,
> WSE2, WSE3, WSE4), but I want to restrict a vector I am going to sample
> from
> to only those columns that are not correspond to GROUP 1 variables I have
> sampled. That is, I want to sample a column from WSE3, WSE4  Columns
> corresponding to CSE1 and CSE2 (i.e., WSE1, WSE2) need to be dropped.
>
> How can I do this? what if I want to repeat this whole process (drawing 2
> random columns from CSE1, CSE2, CSE3, and CSE4 first, AND then another
> random column from WSE1, WSE2, WSE3, and WSE4) for 1000 times. any ideas?
>
> Many thanks in advance!!
>
> --
> View this message in context:
> http://r.789695.n4.nabble.com/Sampling-problem-tp3043804p3043804.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] Anyone can help with this question

[Henrique Dallazuanna]

Try this:

DF[sample(seq(nrow(DF)), 1000),]

Where DF is your data

On Tue, Nov 16, 2010 at 3:13 PM, Mariana  wrote:

>
> Hi there:
> I am a total beginner in R, and I have a simple question:
> I have a table with thousands of lines that represent locations, and two
> columns: latitude and longitude. I need to randomly sample 1000 lines. How
> do I do it? I know the command "sample", but it samples elements
> independently, not lines.
> If there is a better place for me to ask that type of question, please let
> me know.
> Thanks a lot,
> Mariana (from Brazil)
>
> --
> View this message in context:
> http://r.789695.n4.nabble.com/Anyone-can-help-with-this-question-tp3045265p3045265.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.
>



-- 
Henrique Dallazuanna
Curitiba-Paraná-Brasil
25° 25' 40" S 49° 16' 22" O

[[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] help using function merge_all()

[arjun]

Hi Josh,
  Thanks for your reply.

  While I am loading the package, it says package plyr is
required and I tried to install plyr package but I am unable to do so
because it is giving an error message as shown below.

Please help me.




library(reshape)
Loading required package: plyr
Error: package 'plyr' could not be loaded
In addition: Warning message:
In library(pkg, character.only = TRUE, logical.return = TRUE, lib.loc =
lib.loc) :
  there is no package called 'plyr'

> utils:::menuInstallPkgs()
trying URL 'http://cran.mtu.edu/bin/windows/contrib/2.11/plyr_1.1.zip'
Error in download.file(url, destfile, method, mode = "wb", ...) : 
  cannot open URL
'http://cran.mtu.edu/bin/windows/contrib/2.11/plyr_1.1.zip'
In addition: Warning message:
In download.file(url, destfile, method, mode = "wb", ...) :
  cannot open: HTTP status was '404 Not Found'
Warning in download.packages(pkgs, destdir = tmpd, available = available,  :
  download of package 'plyr' failed

-- 
View this message in context: 
http://r.789695.n4.nabble.com/help-using-function-merge-all-tp3045173p3045537.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] Anyone can help with this question

[Jonathan P Daily]

assuming your data takes the form of

locationlatitudelongitude
string  num num
string2 num num

try:

sub <- dat[sample.int(length(dat$location), 1000),]
--
Jonathan P. Daily
Technician - USGS Leetown Science Center
11649 Leetown Road
Kearneysville WV, 25430
(304) 724-4480
"Is the room still a room when its empty? Does the room,
 the thing itself have purpose? Or do we, what's the word... imbue it."
 - Jubal Early, Firefly

r-help-boun...@r-project.org wrote on 11/16/2010 12:13:57 PM:

> [image removed] 
> 
> [R] Anyone can help with this question
> 
> Mariana 
> 
> to:
> 
> r-help
> 
> 11/16/2010 02:22 PM
> 
> Sent by:
> 
> r-help-boun...@r-project.org
> 
> 
> Hi there:
> I am a total beginner in R, and I have a simple question:
> I have a table with thousands of lines that represent locations, and two
> columns: latitude and longitude. I need to randomly sample 1000 lines. 
How
> do I do it? I know the command "sample", but it samples elements
> independently, not lines. 
> If there is a better place for me to ask that type of question, please 
let
> me know.
> Thanks a lot,
> Mariana (from Brazil)
> 
> -- 
> View this message in context: http://r.789695.n4.nabble.com/Anyone-
> can-help-with-this-question-tp3045265p3045265.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] Extending a plot in a loop (With attachment)

[Peter Alspach]

Tena koe Nasrin

Try points() instead of plot() in your second and subsequent calls to plot().  
points() and lines() adds to the current plot by default.  Of course you may 
have difficulties with setting the x and y limits by that's another matter.

HTH 

Peter Alspach

> -Original Message-
> From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-
> project.org] On Behalf Of Nasrin Pak
> Sent: Wednesday, 17 November 2010 7:24 a.m.
> Cc: r-help@r-project.org
> Subject: [R] Extending a plot in a loop (With attachment)
> 
> My problem is that I have a data set for every day of measurement in a
> seperate file and I want to plot one parameter of the data for all the
> days in one graph. I tried to use for loop but only the last data
> remains in the program memory, I don`t know how to plot each day`s data
> continusly after the others(or how to extending the x axis.) Would you
> please help me with it?
> 
> This a plot for one day:
> 
> radiation.data
> <-read.table("C:/updated_CFL_Rad_files/2008/RAD_2008_JD101_0410.dat",
> header = TRUE,sep = ",", quote = " ", dec = ".")
> > attach(radiation.data)
> The following object(s) are masked from 'radiation.data (position 3)':
> 
> Batt_avg, Batt_st, Day, Hour, Kdown_avg, Kdown_st, LW_in, LW_in_st,
> Minute, Month, PanelT_avg, PanelT_st, PAR_avg, PAR_st, Sec,
> Tcase_avg, Tcase_st, Tdome_avg, Tdome_st, Thermopile_avg,
> Thermopile_st, Tuv_avg, Tuv_st, Uva_avg, Uva_st, Uvb_avg, Uvb_st,
> Year
> > names(radiation.data)
>  [1] "Year"   "Month"  "Day""Hour"
>  [5] "Minute" "Sec""Batt_avg"   "PanelT_avg"
>  [9] "Batt_st""PanelT_st"  "Kdown_avg"
> "Thermopile_avg"
> [13] "Tcase_avg"  "Tdome_avg"  "LW_in"  "PAR_avg"
> [17] "Tuv_avg""Uvb_avg""Uva_avg""Kdown_st"
> [21] "Thermopile_st"  "Tcase_st"   "Tdome_st"   "LW_in_st"
> [25] "PAR_st" "Tuv_st" "Uvb_st" "Uva_st"
> 
>  plot(((PAR_avg*0.216)/Uvb_avg),
> main="Par/UVB",xlab="minutes",ylab="Par/UVB")
> 
> 
> and this is the algorithm I tried  for plotting all the data in one
> plot:
> 
> x<- matrix( list.files("C:/updated_CFL_Rad_files",full=TRUE)) # putting
> all data sets in a matrix
>   for(i in 1:100) {
>   if(i < 101) next
>  radiation.data <-read.table(x[i], header = TRUE,sep = ",", quote =
> " ", dec = ".")
>  attach(radiation.data)
>  plot(i*Hour*60+Minute,PAR_avg,main="PAR",xlab="Hour",ylab="Par")
> dev.print(device=postscript, "C:/graph5.eps", onefile=FALSE,
> horizontal=FALSE)
>}
> The plot I see is the last file's plot, I don't know how to keep
> previous data and continue within the same plot.
> 
> * I have attached a sample of the data to this email.
> 
> --
> Sincerely
> 
> Nasrin  Pak
> 
> 
> 
> 
> 
> 
> 
> --
> Sincerely
> 
> Nasrin  Pak

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] help using function merge_all()

[arjun]

Hi Josh,
 I thought of giving up and started writing code in excel using
VBA but then I saw ur message and gave a try in R. I got it, the location
from which I am calling the package "plyr" is not working out so called from
different CRANmirror location and then I was able load the package. now the
merge_all function is working fine. Thanks for your help.  
-- 
View this message in context: 
http://r.789695.n4.nabble.com/help-using-function-merge-all-tp3045173p3045561.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] Sweave: \Sexpr and variables with special chars

[Ralf B]

I am using \Sexpr to include a variable in a title of a Sweave document:

\documentclass[a4paper]{article}
<>=
#mytitlevar <- "Stuff" # case 1, everything is find
mytitlevar <- "Stuff_first" # case 2, f is turned into sub-text
@
\title{MyTitle: \\ \Sexpr{mytitlevar} }
\begin{document}
\maketitle
\end{document}

When doing this, the variable seems to be subject to interpretation by
LaTeX. The variable in case #2 causes the 'f' of 'Stuff_first' to be
printed as sub-text because of the leading underscore. Is there a way
to turn the variable value (the text) into a form so that it is not
interpreted and/or Sweave? I understand that this perhaps more of a
LaTeX question than an R question...

Thanks,
Ralf

__
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 fitting spatial glmm with correlated random effects

[Elijah DePalma]

Greetings,

May you please suggest a package or function to use for fitting a GLMM
(generalized linear mixed model) with spatially correlated random effects?

Thank you,
Elijah DePalma

[[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] Re : interpretation of coefficients in survreg AND obtaining the hazard function for an individual given a set of predictors

[Thomas Lumley]

A coefficient of -0.4 means that survival times are multiplied by
exp(-0.4), that is, people survival only 67% as long.

-thomas


On Wed, Nov 17, 2010 at 4:32 AM, Vincent Vinh-Hung  wrote:
> Thanks for sharing the questions and responses!
>
> Is it possible to appreciate how much the coefficients matter in one
> or the other model?
> Say, using Biau's example, using coxph, as.factor(grade2 ==
> "high")TRUE gives hazard ratio 1.27 (rounded).
> As clinician I can grasp this HR as 27% relative increase. I can
> relate with other published results.
> With survreg the Weibull model gives a coefficient -0.4035245: is it
> feasible or meaningful to translate it to HR?
>
> Thanks in advance,
>
> Vincent Vinh-Hung
> Radiation Oncology,
> Geneva University Hospitals
>
> On Sun, Nov 14, 2010 at 6:51 AM, Biau David  wrote:
>> Dear R help list,
>>
>> I am modeling some survival data with coxph and survreg (dist='weibull') 
>> using
>> package survival. I have 2 problems:
>>
>> 1) I do not understand how to interpret the regression coefficients in the
>> survreg output and it is not clear, for me, from ?survreg.objects how to.
>>
>> Here is an example of the codes that points out my problem:
>> - data is stc1
>> - the factor is dichotomous with 'low' and 'high' categories
>>
>> slr <- Surv(stc1$ti_lr, stc1$ev_lr==1)
>>
>> mca <- coxph(slr~as.factor(grade2=='high'), data=stc1)
>> mcb <- coxph(slr~as.factor(grade2), data=stc1)
>> mwa <- survreg(slr~as.factor(grade2=='high'), data=stc1, dist='weibull',
>> scale=0)
>> mwb <- survreg(slr~as.factor(grade2), data=stc1, dist='weibull', scale=0)
>>
>>> summary(mca)$coef
>>                                                             coef
>> exp(coef)      se(coef)         z                      Pr(>|z|)
>> as.factor(grade2 == "high")TRUE 0.2416562  1.273356     0.2456232
>> 0.9838494      0.3251896
>>
>>> summary(mcb)$coef
>>                                       coef             exp(coef)
>> se(coef)             z                     Pr(>|z|)
>> as.factor(grade2)low -0.2416562 0.7853261     0.2456232     -0.9838494
>> 0.3251896
>>
>>> summary(mwa)$coef
>> (Intercept)     as.factor(grade2 == "high")TRUE
>> 7.9068380       -0.4035245
>>
>>> summary(mwb)$coef
>> (Intercept)     as.factor(grade2)low
>> 7.5033135       0.4035245
>>
>>
>> No problem with the interpretation of the coefs in the cox model. However, i 
>> do
>> not understand why
>> a) the coefficients in the survreg model are the opposite (negative when the
>> other is positive) of what I have in the cox model? are these not the log(HR)
>> given the categories of these variable?
>
> No. survreg() fits accelerated failure models, not proportional
> hazards models.   The coefficients are logarithms of ratios of
> survival times, so a positive coefficient means longer survival.
>
>
>> b) how come the intercept coefficient changes (the scale parameter does not
>> change)?
>
> Because you have reversed the order of the factor levels.  The
> coefficient of that variable changes sign and the intercept changes to
> compensate.
>
>
>> 2) My second question relates to the first.
>> a) given a model from survreg, say mwa above, how should i do to extract the
>> base hazard and the hazard of each patient given a set of predictors? With 
>> the
>> hazard function for the ith individual in the study given by  h_i(t) =
>> exp(\beta'x_i)*\lambda*\gamma*t^{\gamma-1}, it doesn't look like to me that
>> predict(mwa, type='linear') is \beta'x_i.
>
> No, it's beta'x_i for the accelerated failure parametrization of the
> Weibull.  In terms of the CDF
>
> F_i(t) = F_0( exp((t+beta'x_i)/scale) )
>
> So you need to multiply by the scale parameter and change sign to get
> the log hazard ratios.
>
>
>> b) since I need the coefficient intercept from the model to obtain the scale
>> parameter  to obtain the base hazard function as defined in Collett
>> (h_0(t)=\lambda*\gamma*t^{\gamma-1}), I am concerned that this coefficient
>> intercept changes depending on the reference level of the factor entered in 
>> the
>> model. The change is very important when I have more than one predictor in 
>> the
>> model.
>
> As Terry Therneau pointed out recently in the context of the Cox
> model, there is no such thing as "the" baseline hazard.  The baseline
> hazard is the hazard when all your covariates are equal to zero, and
> this depends on how you parametrize.  In mwa, zero is grade2="low", in
> mwb, zero is grade2="high", so the hazard at zero has to be different
> in the two cases.
>
>     -thomas
>
> --
> Thomas Lumley
> Professor of Biostatistics
> University of Auckland
>
> __
> 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.
>



-- 
Thomas Lumley
Professor of Biostatistics
University of Auckland

_

Re: [R] Computing Rolling Average

[Ray Brownrigg]

On Wed, 17 Nov 2010, Paolo Rossi wrote:
> Hi,
> Can anyone suggest a  clever way to compute a rolling weekly average of the
> columns in a matrix?  The column bit is straightforward – use apply given a
> function which does what you want on a column. With regard to a particular
> column, the obvious way is to run a for loop indexing the last 7 days and
> computing the average . I simply would like to know if there is a  better /
> quicker way.
>
>
>
> Code:
> Given a,
>
> > a= array(1:100, dim = c(17,2))
> > a
>
>   [,1] [,2]
>  [1,]1   18
>  [2,]2   19
>  [3,]3   20
>  [4,]4   21
>  [5,]5   22
>  [6,]6   23
>  [7,]7   24
>  [8,]8   25
>  [9,]9   26
> [10,]   10   27
> [11,]   11   28
> [12,]   12   29
> [13,]   13   30
> [14,]   14   31
> [15,]   15   32
> [16,]   16   33
> [17,]   17   34
> one needs to start computing the average from obs 7 s (at obs 7 you have a
> full week to compute the average) and compute the rolling weekly average
> from day 7 onwards
> Results will look like b
>   [,1] [,2]
>  [1,]4   14
>  [2,]5   21
>  [3,]6   22
>  [4,]7   23
>  [5,]8   24
>  [6,]9   25
>  [7,]   10   26
>  [8,]   11   27
>  [9,]   12   28
> [10,]   13   29
> Thanks in advance,
>
I don't see how an average of 7 numbers all 18 or greater can be 14, as in your
result[1, 2], unless you have mis-stated the question.

Anyway, try:
apply(a, 2, function(x) {cx <- cumsum(x); N <- length(x); (cx[7:N] - c(0, 
cx[1:(N-7)]))/7}

HTH
Ray Brownrigg

> Paolo

__
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] Offset in glm poisson using R vs Exposure in Stata

[Columbine Caroline Waring]

Ben, 
Thank you, it WAS a typo of sorts.

Officially I tried:
> glm(count~md+ms+rf+sg+offset(log(Eff)),family=poisson,data=DepthHabGen)
> glm(count~md+ms+rf+sg, offset=(log(Eff)),family=poisson,data=DepthHabGen)
(which of course are the same as eachother)

> glm(count~md+ms+rf+sg, offset=(Eff),family=poisson,data=DepthHabGen)
> glm(count~md+ms+rf+sg+offset(Eff),family=poisson,data=DepthHabGen)
(which are also the same between themselves, yet wrong compared to the STATA 
model)

Additionally, given the text you found on stata website, which I am familiar 
with, I also tried:

> glm(count~md+ms+rf+sg, offset=(exp(Eff)),family=poisson,data=DepthHabGen)

> glm(count~md+ms+rf+sg+offset(exp(Eff)),family=poisson,data=DepthHabGen)
(which still might be the solution however R issues the following response:
Error: no valid set of coefficients has been found: please supply starting 
value)

Which, according to other help postings might be able to be forced to function. 
 And I quote from 2008 message 
[R] glm error message when using family Gamma(link="inverse")"You might be able 
to calm glm()'s frayed nerves by supplying some decent 
start values as it is asking. Possibly starting with a subset of 
variables, using the coef()s for the subset and setting the others to zero 
when you attempt to fit all together will be enough to push it past its 
sticking point."

I'm not sure how to go about employing this:
Using, coef() will bring up a list of all coeff, which at this point having 
taken the inverse log or exponential of the variable "Eff", are all very large!

(Eff: is a set of real numbers representing the mean of effort in days 1-209, 
for each count 
 

Also, combined with other postings  and I gather it may be possible to use the 
start() or etastart() functions 

>From [R] documentation, Fitting Generalized linear Models

start

starting values for the parameters in the linear predictor.
etastart

starting values for the linear predictor.
Which I don't understand how to do?

Or perhaps change the link? What would be recommended?

Thank you.

> To: r-h...@stat.math.ethz.ch
> From: bbol...@gmail.com
> Date: Tue, 16 Nov 2010 13:16:20 +
> Subject: Re: [R] Offset in glm poisson using R vs Exposure in Stata
> 
> Columbine Caroline Waring  hotmail.com> writes:
> 
> > I am hoping to find someone who uses both R and program Stata for GLMs.
> [snip]
> > What I have is the code from Stata  and am trying to reproduce the same
> analysis in R - my program of choice.
> > 
> > . glm count md ms rf sg, family(poisson)
> > exposure(effort) eform
> > 
> > I am lost at the point of finding the equivalent code for 'exposure'. 
> > 
> > Having looked at a few forums and 'googled'. I thought 'offset', used as   
> offset=(log(Eff))   or the
> > equivalent+offset(log(Eff))   would produce the desired effect.
> > 
> > Incidentally my code was:  
>   glm(Count~md+ms+rf+sg+offset(Eff),family=poisson,data=DepthHabGen)
> > 
> > (Making use of glm{stats})
> > 
> 
>   Based on your discussion above did you mean
> 
> glm(count~md+ms+rf+sg+offset(log(Eff)),
>family=poisson,data=DepthHabGen)
> 
>  ? is this just a typo ?
> 
> according to http://www.stata.com/help.cgi?glm
> 
>  exposure(varname) include ln(varname) in model with coefficient
>   constrained to 1
>  offset(varname)   include varname in model with coefficient
>   constrained to 1
> 
> 
> > However, offset does not seem to be equivalent to 'exposure' in Stata. As
> coefficients and log likelhood
> > estimates differ.
> > 
> > So I asked the following questions:
> > 
> > 1. Do both programs produce the same results without
> > 'exposure' i.e. glm models 
> > 
> > Yes,  log likelihoods and coefficients are the same.
> > 
> > 2. How about using the unintuitive non logged " offset=Eff" ?
> > 
> > Coefficients and log likelihoods still differ.
> 
>   You're not doing anything obviously wrong.  Are you sure your
> "effort" and "Eff" variables are the same, i.e. nothing got mangled
> moving to R?
> 
>   I don't use Stata, perhaps someone else can try.  Posting a small
> reproducible example would be helpful.
> 
> __
> 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] format secondary axis for dates

[David Winsemius]

On Nov 16, 2010, at 2:32 PM, Jannis wrote:


Dear List,


this may be a Newbi question and may have been asked several times,  
but i am too stupid to find the posts.



I have a plot of values against POSIXct time steps. If I want to add  
a second x axis to the top margin of the plot, only numbers are at  
the tickmarks. Is there a straightforward way to specify the format  
to convert them to character representations (similar to the labels  
at the bottom)?


x<-as.POSIXct(1:1000*(60^2),origin='01-01-1970')


The origin argument was incorrect and you got the wrong starting point:
> x[1]
[1] "0001-01-19 01:00:00 EST"



y=rnorm(1000)
plot(x,y)
axis(3)


This gives ticks at hourly intervals:

axis(3, labels=format(x, "%Y-%m-%d"), at=x)

Whereas this is probably what you want (after correcting the origin to  
"1970-01-01"


axis.POSIXct(3, x=x, labels=TRUE )

--
David.






Thanks for your help
Jannis





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


David Winsemius, MD
West Hartford, CT

__
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] Computing Rolling Average

[David Winsemius]

On Nov 16, 2010, at 2:33 PM, Paolo Rossi wrote:


Hi,
Can anyone suggest a  clever way to compute a rolling weekly average  
of the
columns in a matrix?  The column bit is straightforward – use apply  
given a
function which does what you want on a column. With regard to a  
particular
column, the obvious way is to run a for loop indexing the last 7  
days and
computing the average . I simply would like to know if there is a   
better /

quicker way.


There is a rollmean function in package zoo. See below.





Code:
Given a,

a= array(1:100, dim = c(17,2))
a

 [,1] [,2]
[1,]1   18
[2,]2   19
[3,]3   20
[4,]4   21
[5,]5   22
[6,]6   23
[7,]7   24
[8,]8   25
[9,]9   26
[10,]   10   27
[11,]   11   28
[12,]   12   29
[13,]   13   30
[14,]   14   31
[15,]   15   32
[16,]   16   33
[17,]   17   34
one needs to start computing the average from obs 7 s (at obs 7 you  
have a
full week to compute the average) and compute the rolling weekly  
average

from day 7 onwards
Results will look like b
 [,1] [,2]
[1,]4   14
[2,]5   21


I think you got your first ones out of sequence:

> apply(a, 2, rollmean, k=7)
  [,1] [,2]
 [1,]4   21
 [2,]5   22
 [3,]6   23
 [4,]7   24
 [5,]8   25
 [6,]9   26
 [7,]   10   27
 [8,]   11   28
 [9,]   12   29
[10,]   13   30
[11,]   14   31


--
David.

[3,]6   22
[4,]7   23
[5,]8   24
[6,]9   25
[7,]   10   26
[8,]   11   27
[9,]   12   28
[10,]   13   29
Thanks in advance,

Paolo



On 15 November 2010 20:59, wangwallace  wrote:



Hey,

I am hoping someone can help me with a sampling question.

I have a data frame of 8 variables (the first column is the  
subjects' id):


  SubIDCSE1 CSE2 CSE3 CSE4 WSE1 WSE2 WSE3 WSE4
1  6  5   6   2  6  2 
2   4
2  6  4   7   2  6  6 
2   3
3  5  5   5   5  5  5 
4   5
4  5  4   3   4  4  4 
5   2
5  5  6   7   5  6  4 
4   1
6  5  4   3   6  4  3 
7   3
7  3  6   6   3  6  5 
2   1
8  3  6   6   3  6  5 
4   7


the 6 variables are categorized into two groups with CSE1, CSE2,  
CSE3, and

CSE4 in one group and the rest in another group.


sample(data[,2:4],2,replace=FALSE)


 CSE1 CSE2
1  65
2  64
3  55
4  54
5  56
6  54
7  36
8  36

Now I want to sample 1 column from another group of variables  
(i.e., WSE1,
WSE2, WSE3, WSE4), but I want to restrict a vector I am going to  
sample

from
to only those columns that are not correspond to GROUP 1 variables  
I have

sampled. That is, I want to sample a column from WSE3, WSE4  Columns
corresponding to CSE1 and CSE2 (i.e., WSE1, WSE2) need to be dropped.

How can I do this? what if I want to repeat this whole process  
(drawing 2
random columns from CSE1, CSE2, CSE3, and CSE4 first, AND then  
another
random column from WSE1, WSE2, WSE3, and WSE4) for 1000 times. any  
ideas?


Many thanks in advance!!

--
View this message in context:
http://r.789695.n4.nabble.com/Sampling-problem-tp3043804p3043804.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.


David Winsemius, MD
West Hartford, CT

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