[R] plot xaxp issue

2013-01-07 Thread Elaine Kuo
Hello,

I have data of Body length and Body weight of 8 boys and 7 girls.

I want to draw the plot of Body length (for X) and Body weight (for Y)
based on sex.
Then the two plots want to be overlapped for comparison.

I used the code below but found the unit length of X axis of boy and
girl plot are not the same.
For instance, the length between 0 and 1 of boy plot is larger than
that in girl plot.
The same thing happened to Y axis as well.
(In other words, though axap and yaxp were set to be the same, the
display were not the same.)

Please kindly advise correction of the code.
Thank you.

Elaine

# plot code
boy<-read.csv("H:/boy_data.csv",header=T)
girl<-read.csv("H:/girl_data.csv",header=T)
par(mai=c(1.03,1.03,0.4,0.4))

plot(boy$body_length, boy$body_weight,
xlab=" body_length (cm)",
ylab=" body_weight  ( kg )",
xaxp=c(0,200,4),
yaxp=c(0,100,4),
type="p",
pch=1,lwd=1.0,
cex.lab=1.4, cex.axis=1.2,
font.axis=2,
cex=1.5,
las=1,
bty="l",col="firebrick3")

boyline<-lm(body_weight ~ body_length, boy)
summary(boyline)
abline(boyline,col="firebrick3",lwd=2)

#~~~
  # graph
par(mai=c(1.03,1.03,0.4,0.4))

par(new=T)

plot(girl$body_length, girl$body_weight,
xlab=" body_length (cm)",
ylab=" body_weight  ( kg )",
xaxp=c(0,200,4),
yaxp=c(0,100,4),
type="p",
pch=1,lwd=1.0,
cex.lab=1.4, cex.axis=1.2,
font.axis=2,
cex=1.5,
las=1,
bty="l",col="saddlebrown")


girlline<-lm(body_weight~ body_length, girl)
summary(girlline)
abline(girlline,col="saddlebrown",lwd=2)

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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 xaxp issue

2013-01-07 Thread PIKAL Petr
Hi

Instead of two plots with par(new = TRUE) try to put boys and girls together 
(quite natural thing, they will be pleased 8-)

together <- rbind(boy, girl)
together$sex <- factor(rep(c("boy", "girl"), c(8,7)))

plot(together$body_length, together$body_weight, , 
col=c("firebrick3","saddlebrown")[as.numeric(together$sex)], )

Regards
Petr


> -Original Message-
> From: r-help-boun...@r-project.org [mailto:r-help-bounces@r-
> project.org] On Behalf Of Elaine Kuo
> Sent: Monday, January 07, 2013 9:27 AM
> To: r-help@r-project.org
> Subject: [R] plot xaxp issue
> 
> Hello,
> 
> I have data of Body length and Body weight of 8 boys and 7 girls.
> 
> I want to draw the plot of Body length (for X) and Body weight (for Y)
> based on sex.
> Then the two plots want to be overlapped for comparison.
> 
> I used the code below but found the unit length of X axis of boy and
> girl plot are not the same.
> For instance, the length between 0 and 1 of boy plot is larger than
> that in girl plot.
> The same thing happened to Y axis as well.
> (In other words, though axap and yaxp were set to be the same, the
> display were not the same.)
> 
> Please kindly advise correction of the code.
> Thank you.
> 
> Elaine
> 
> # plot code
> boy<-read.csv("H:/boy_data.csv",header=T)
> girl<-read.csv("H:/girl_data.csv",header=T)
> par(mai=c(1.03,1.03,0.4,0.4))
> 
> plot(boy$body_length, boy$body_weight,
> xlab=" body_length (cm)",
> ylab=" body_weight  ( kg )",
> xaxp=c(0,200,4),
> yaxp=c(0,100,4),
> type="p",
> pch=1,lwd=1.0,
> cex.lab=1.4, cex.axis=1.2,
> font.axis=2,
> cex=1.5,
> las=1,
> bty="l",col="firebrick3")
> 
> boyline<-lm(body_weight ~ body_length, boy)
> summary(boyline)
> abline(boyline,col="firebrick3",lwd=2)
> 
> #~~~
>   # graph
> par(mai=c(1.03,1.03,0.4,0.4))
> 
> par(new=T)
> 
> plot(girl$body_length, girl$body_weight,
> xlab=" body_length (cm)",
> ylab=" body_weight  ( kg )",
> xaxp=c(0,200,4),
> yaxp=c(0,100,4),
> type="p",
> pch=1,lwd=1.0,
> cex.lab=1.4, cex.axis=1.2,
> font.axis=2,
> cex=1.5,
> las=1,
> bty="l",col="saddlebrown")
> 
> 
> girlline<-lm(body_weight~ body_length, girl)
> summary(girlline)
> abline(girlline,col="saddlebrown",lwd=2)
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-
> guide.html
> and provide commented, minimal, self-contained, reproducible code.

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


Re: [R] how to aggregate T-test result in an elegant way?

2013-01-07 Thread Yao He
Hi, arun
I'm so sorry for that isn't helpful.
One of question is that I don't know how  to subset a small part as it
is a 3-dimension array so I just show the structure of that.
 I tried  dput()  to a file , then what should I do for subsetting it?

Another question is :
My rawdata is a "melt" dataframe like that:
IID O2  variablevalue
1   TWF2H5  13% EW.INCU 49.38
2   TWF2H6  13% EW.INCU 48.02
3   TWF2H19 13%  EW.INCU51.44
280 TWF2H10113% EW.17.5 42.26
281 TWF2H10513%  EW.17.5 43.52
282 TWF2H10613% EW.17.5 42.83
472 TWF2N10221% EW.17.5 45.97
473 TWF2N10421%  EW.17.543.32
474 TWF2N10621% EW.17.5 48.63
689 TWF2N2  21%  EMW19.57
690 TWF2N6  21%  EMW18.07
691 TWF2N10 21% EMW 15.4
491 TWF2H5   13%EMW 15.61
492 TWF2H6  13% EMW 13.41
493 TWF2H19 13% EMW 14.03
199 TWF2N2  21% EW.INCU 48.69
200 TWF2N6  21% EW.INCU 50.52
201 TWF2N10 21% EW.INCU 42.04

if you meet a t-test task as I described  , is that generate a
high-dimension array  a good way ?
Thank you!

Yao He
2013/1/7 arun :
> HI,
> I tried to create an example dataset (as you didn't provide the data).
> set.seed(25)
> a<-array(sample(1:50,60,replace=TRUE),dim=c(2,10,3))
> dimnames(a)[[1]]<-c("13%","21%")
> dimnames(a)[[2]]<-paste("TWF2H",101:110,sep="")
> dimnames(a)[[3]]<-c("EW.INCU","EW.17.5","EMW")
>
>
> str(a)
> # int [1:2, 1:10, 1:3] 21 35 8 45 7 50 32 17 4 15 ...
>  #- attr(*, "dimnames")=List of 3
>   #..$ : chr [1:2] "13%" "21%"
>   .#.$ : chr [1:10] "TWF2H101" "TWF2H102" "TWF2H103" "TWF2H104" ...
>   #..$ : chr [1:3] "EW.INCU" "EW.17.5" "EMW"
>
> res<-lapply(lapply(seq_len(dim(a)[3]),function(i) 
> t.test(a[dimnames(a)[[1]][1],,i],a[dimnames(a)[[1]][2],,i])),function(x) 
> data.frame(mean=x$estimate,p.value=x$p.value))
> res1<-do.call(rbind,res)
>   row.names(res1)[grep("mean of 
> x",row.names(res1))]<-gsub("(.*\\.).*$","\\113%",row.names(res1)[grep("mean 
> of x",row.names(res1))])
>  row.names(res1)[grep("mean of 
> y",row.names(res1))]<-gsub("(.*\\.).*$","\\121%",row.names(res1)[grep("mean 
> of y",row.names(res1))])
> res1
> #mean   p.value
> #EW.INCU.13% 22.3 0.2754842
> #EW.INCU.21% 29.3 0.2754842
> #EW.17.5.13% 20.5 0.4705772
> #EW.17.5.21% 16.0 0.4705772
> #EMW.13% 23.9 0.9638679
> #EMW.21% 24.2 0.9638679
> A.K.
>
>
>
>
> - Original Message -
> From: Yao He 
> To: arun 
> Cc: R help 
> Sent: Sunday, January 6, 2013 11:21 PM
> Subject: Re: [R] how to aggregate T-test result in an elegant way?
>
> Thank you,it is really helpful everytime.
>
> I didn't provide any example data because I thought it is just a
> question of how to report t.test() result in R.
> However,as you say,it is better to show more details for finding an elegant 
> way
>
> In fact  I generate a 3-dimension array like that:
> str(a)
> num [1:2, 1:245, 1:3] 47.5 NA 48.9 NA 47.5 ...
> - attr(*, "dimnames")=List of 3
>   ..$ : chr [1:2] "13%" "21%"
>   ..$ : chr [1:245] "TWF2H101" "TWF2H105" "TWF2H106" "TWF2H110" ...
>   ..$ : chr [1:3] "EW.INCU" "EW.17.5" "EMW"
>
> I want to do two sample mean t-test between 13% and 21% for each
> variable "EW.INCU" "EW.17.5" "EMW".
>
> So I try these codes:
> variable<-dimnames(a)[[3]]
>   O2<-dimnames(a)[[1]]
>   for (i in variable) {
> print(i)
> print(O2[1])
> print(O2[2])
> print(t.test(a[O2[1],,i],a[O2[2],,i],na.rm=T))
> }
>
> I don't think it is an elegant way and I am inexperience to report raw result.
> Could you give me more help?
>
> Yao He
>
> 2013/1/7 arun :
>> Hi,
>> You didn't provide any example data.  So, I am not sure whether this helps.
>>
>> set.seed(15)
>> dat1<-data.frame(A=sample(10:20,5,replace=TRUE),B=sample(18:28,5,replace=TRUE),C=sample(25:35,5,replace=TRUE),D=sample(20:30,5,replace=TRUE))
>>  res<-lapply(lapply(seq_len(ncol(dat2)),function(i) 
>> t.test(dat2[,i],dat1[,1],paired=TRUE)),function(x) 
>> data.frame(meanDiff=x$estimate,p.value=x$p.value))# paired
>> names(res)<-paste("A",LETTERS[2:4],sep="")
>> res<- do.call(rbind,res)
>> res
>>   # meanDiff p.value
>> #AB  9.4 0.021389577
>> #AC 15.0 0.002570261
>> #AD 10.6 0.003971604
>>
>>
>> #or
>> res1<-lapply(lapply(seq_len(ncol(dat2)),function(i) 
>> t.test(dat2[,i],dat1[,1],paired=FALSE)),function(x) 
>> data.frame(mean=x$estimate,p.value=x$p.value))
>> names(res1)<-paste("A",LETTERS[2:4],sep="")
>> res1<-do.call(rbind,res1)
>> row.names(res1)[grep("mean of 
>> y",row.names(res1))]<-gsub("(.*\\.).*","\\1A",row.names(res1)[grep("mean of 
>> y",row.names(res1))])
>> row.names(res1)[grep("mean of 
>> x",row.names(res1))]<-gsub("(\\w)(\\w)(\\.).*","\\1\\2\\3\\2",row.names(res1)[grep("mean
>>  of x",row.names(res1))])
>> res1
>> # mean  p.value
>> #AB.B 25.2 1.299192e-03
>> #AB.A 15.8 1.299192e-03
>> #AC.C 30.8 5.145519e-05
>> #

[R] round 2.3 to 2.5 and 2.2 to 2.0

2013-01-07 Thread Mat
Hello together,

i want to round some numbers in my data.frame.
How can i round a number to values like 0.5; 1.0; 1.5; 2.0; etc.

It should look like this one

before
2.2 ; 2.3; 2.26; 1.11

after
2.0; 2.5; 2.5; 1.0

thanks.

Mat



--
View this message in context: 
http://r.789695.n4.nabble.com/round-2-3-to-2-5-and-2-2-to-2-0-tp4654813.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] plot x-axis DateTime NOT evenly spaced

2013-01-07 Thread ishi soichi
R-64 latest

Hi. I am trying to plot a set of csv data, which looks like

> head(interval)
   date inteval
1 2012-07-01 00:57:54 +0900 156
2 2012-07-01 01:07:41 +0900 587
3 2012-07-01 01:09:31 +0900 110
4 2012-07-01 01:18:42 +0900 551
5 2012-07-01 01:39:01 +09001219
6 2012-07-01 01:40:40 +0900  99

as you can see, more than one event happens each day, and they are not
evenly spaced.  Obviously hours, minutes and seconds are important for the
plot.

I tried

interval$date <- as.Date(interval$date, "%Y-%m-%d %H:%M:%S +0900")

but this chops the time off.

Could anyone show me how to plot data with x values as Date(or Time)
objects?

soichi

[[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] package metafor: error when setting 'col' and 'at' for a forest plot

2013-01-07 Thread Viechtbauer Wolfgang (STAT)
Dear Brian,

At the moment, the various forest() functions are not meant to accept 'col' as 
an argument. While it is indeed possible to specify a 'col' argument, it will 
be passed on via the ... argument to further functions within forest() and this 
is where things can go awry. 

To give a more technical explanation, let's look at your example. First, note 
that the highest 'at' value you specified is log(10), which is approximately 
2.30. However, the highest upper CI bound in your dataset is 2.9. Therefore, 
the forest() function wants to draw a right-pointed arrow for that study using 
the polygon() function. The code within forest() looks something like this:

polygon(, col="black", ...)

When you specify a 'col' argument, it is passed on to polygon() via ... -- but 
col="black" is already hard-coded in forest() and so you get the "formal 
argument "col" matched by multiple actual arguments" error.

The issue of changing colors in the forest() functions has come up before:

https://stat.ethz.ch/pipermail/r-help/2011-August/287788.html

In fact, the problems associated with allowing the user to specify colors in 
high-level plotting functions, such as forest(), was nicely discussed in a 
recent article by Paul Murrell in the R Journal:

http://journal.r-project.org/archive/2012-2/RJournal_2012-2_Murrell.pdf

I may add functionality to the forest() functions to handle a 'col' argument 
for specifying row colors, but I will have to evaulate carefully how well that 
will work. Also, I would then need to change the current 'col' argument in the 
forest.rma() function -- breaking backwards compatability :/

Back to your example. You can actually get it to work if you use:

forest(forest$OR, ci.lb=forest$Low, ci.ub=forest$High,  at=log(c(.05, .25, 1, 
20)), slab=forest$SNP, atransf=exp)

in which case the values of 'at' encompass the CI bounds of all studies. Then 
the forest() function does not use polygon(), but segments() and here 'col' is 
not hard-coded. One (probably) unintended consequence of using 'col' this way 
is that the x axis label is also given a different color. 

A few other notes:

1) Calling your data frame 'forest' is probably not a good idea -- since that 
is the name of a function.

2) Your data appear to be (raw -- not log transformed) odds ratios and the 
corresponding CI bounds. So, you should not be using the atransf arguments -- 
or you would be transforming your odds ratios to exp(ORs). Probably closer to 
what you want is:

dat <- read.table(header=TRUE, text="
SNP   Group   High   Low   OR
rs1137101   A   1.21   0.87   1.03
rs1137101   B   2.11   1.21   1.6
rs1137101   C   2.9   1.42   2.03
rs1042522   A   1.12   0.84   0.97
rs1042522   B   1.15   0.79   0.95
rs1042522   C   0.92   0.5   0.7
rs1625895   A   1.14   0.76   0.93
rs1625895   B   1.15   0.75   0.93
rs1625895   C   NA   NA   NA
ACEI/D   A   1.55   0.79   1.11
ACEI/D   B   1.25   0.76   0.98
ACEI/D   C   0.85   0.41   0.59
")

forest(dat$OR, ci.lb=dat$Low, ci.ub=dat$High, col=c(1,2,3), at=c(.25, 1, 2, 3), 
slab=dat$SNP, refline=1, xlab=" ", xlim=c(-1.5,5.5))
mtext("Odds Ratio", side=1, line=2.5, adj=0.45)

where I manually add the x axis label, so that it is in black. Again, though, 
using 'col' this way is technically not intended (even though it does work 
here).

I hope this helps!

Best,
Wolfgang

--   
Wolfgang Viechtbauer, Ph.D., Statistician   
Department of Psychiatry and Psychology   
School for Mental Health and Neuroscience   
Faculty of Health, Medicine, and Life Sciences   
Maastricht University, P.O. Box 616 (VIJV1)   
6200 MD Maastricht, The Netherlands   
+31 (43) 388-4170 | http://www.wvbauer.com   

> -Original Message-
> From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org]
> On Behalf Of Brian Z Ring
> Sent: Sunday, January 06, 2013 05:48
> To: 'David Winsemius'; r-help@r-project.org
> Subject: Re: [R] package metafor: error when setting 'col' and 'at' for a
> forest plot
> 
> Good suggestion, and yes, the given examples work. Using a modification of
> one of the supplied examples, where in the final line I only added the
> parameter col=c(1,2,3):
> 
> data(dat.bcg)
> dat <- escalc(measure="RR", ai=tpos, bi=tneg, ci=cpos, di=cneg,
> data=dat.bcg, append=TRUE)
> forest(dat$yi, dat$vi, slab=paste(dat$author, dat$year, sep=", "),
> atransf=exp, at=log(c(.05,.25,1,4,20)), col=c(1,2,3))
> 
> I get a plot similar to what I want to achieve. My code:
> 
> forest(forest$OR, ci.lb=forest$Low, ci.ub=forest$High,  col=c(1,2,3),
> at=log(c(.05, .25, 1, 10)), slab=forest$SNP, atransf=exp)
> 
> differs in the data and using confidence intervals instead of vi (sample
> variance). The parameters being changed in my line of code are fairly
> simple, just color of each row and the x axis tick marks. Setting each of
> these independently with my data works fine, it's only when I try to do
> both
> that it fails. I'm still stuck.
> 
> Here is some example data that should reproduce the error

Re: [R] round 2.3 to 2.5 and 2.2 to 2.0

2013-01-07 Thread Uwe Ligges



On 07.01.2013 09:59, Mat wrote:

Hello together,

i want to round some numbers in my data.frame.
How can i round a number to values like 0.5; 1.0; 1.5; 2.0; etc.

It should look like this one

before
2.2 ; 2.3; 2.26; 1.11

after
2.0; 2.5; 2.5; 1.0



round(2*x)/2

Uwe Ligges


thanks.

Mat



--
View this message in context: 
http://r.789695.n4.nabble.com/round-2-3-to-2-5-and-2-2-to-2-0-tp4654813.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] plot x-axis DateTime NOT evenly spaced

2013-01-07 Thread Uwe Ligges



On 07.01.2013 09:55, ishi soichi wrote:

R-64 latest

Hi. I am trying to plot a set of csv data, which looks like


head(interval)

date inteval
1 2012-07-01 00:57:54 +0900 156
2 2012-07-01 01:07:41 +0900 587
3 2012-07-01 01:09:31 +0900 110
4 2012-07-01 01:18:42 +0900 551
5 2012-07-01 01:39:01 +09001219
6 2012-07-01 01:40:40 +0900  99

as you can see, more than one event happens each day, and they are not
evenly spaced.  Obviously hours, minutes and seconds are important for the
plot.

I tried

interval$date <- as.Date(interval$date, "%Y-%m-%d %H:%M:%S +0900")

but this chops the time off.

Could anyone show me how to plot data with x values as Date(or Time)
objects?



See ?strptime

Best,
Uwe Ligges







soichi

[[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] Changing mtext direction, or using text for the margin?

2013-01-07 Thread Uwe Ligges



On 07.01.2013 07:00, Michael Rennie wrote:


Hi all,

I have read through the archives, but can't find a solution to this
problem.

I need the text direction on "dependent B", plotted in margin 4, to go
top to bottom (opposite what it is now). Here's some sample code:

#plot with mtext example

par(mgp = c(2,1,0), mfrow=c(2,2), las=1, mar=c(2,2,2,2), omi=
c(0.5,0.2,0,0.2))

a<-1:10
b<-7:16
c<-21:30

plot(a~b, ylab="", xlab="", xaxt="n")
axis(1, at=b, labels=FALSE, tick=TRUE)
plot(c~b, ylab="", xlab="", xaxt="n", yaxt="n")
axis(1, at=b, labels=FALSE, tick=TRUE)
axis(4, at=c(22,24,26,28,30), labels=TRUE, tick=TRUE)
plot(a~b, ylab="")
plot(c~b, ylab="", yaxt="n")
axis(4, at=c(22,24,26,28,30), labels=TRUE, tick=TRUE)

mtext(side=1, "independent", outer=TRUE, line=1, padj=1)
mtext(side=2, "dependent A", las=0, outer=TRUE, line=0.25)
mtext(side=4, "dependent B", las=0, outer=TRUE, line=0.25)

I have seen an example in help where i can use text and (srt) to
manipulate this with a single panel plot, but not for a multi-panel
example. Any help would be greatly appreciated.



You cannot do that with mtext, you rather need a dirty hack with text(), 
I believe.


Best,
Uwe Ligges




Cheers,

Mike

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] plot xaxp issue

2013-01-07 Thread Elaine Kuo
Thanks a lot.
Please kindly indicate the meaning of the c(8,7).

> together <- rbind(boy, girl)
> together$sex <- factor(rep(c("boy", "girl"), c(8,7)))

Elaine

On Mon, Jan 7, 2013 at 4:55 PM, PIKAL Petr  wrote:
> Hi
>
> Instead of two plots with par(new = TRUE) try to put boys and girls together 
> (quite natural thing, they will be pleased 8-)
>
> together <- rbind(boy, girl)
> together$sex <- factor(rep(c("boy", "girl"), c(8,7)))
>
> plot(together$body_length, together$body_weight, , 
> col=c("firebrick3","saddlebrown")[as.numeric(together$sex)], )
>
> Regards
> Petr
>
>
>> -Original Message-
>> From: r-help-boun...@r-project.org [mailto:r-help-bounces@r-
>> project.org] On Behalf Of Elaine Kuo
>> Sent: Monday, January 07, 2013 9:27 AM
>> To: r-help@r-project.org
>> Subject: [R] plot xaxp issue
>>
>> Hello,
>>
>> I have data of Body length and Body weight of 8 boys and 7 girls.
>>
>> I want to draw the plot of Body length (for X) and Body weight (for Y)
>> based on sex.
>> Then the two plots want to be overlapped for comparison.
>>
>> I used the code below but found the unit length of X axis of boy and
>> girl plot are not the same.
>> For instance, the length between 0 and 1 of boy plot is larger than
>> that in girl plot.
>> The same thing happened to Y axis as well.
>> (In other words, though axap and yaxp were set to be the same, the
>> display were not the same.)
>>
>> Please kindly advise correction of the code.
>> Thank you.
>>
>> Elaine
>>
>> # plot code
>> boy<-read.csv("H:/boy_data.csv",header=T)
>> girl<-read.csv("H:/girl_data.csv",header=T)
>> par(mai=c(1.03,1.03,0.4,0.4))
>>
>> plot(boy$body_length, boy$body_weight,
>> xlab=" body_length (cm)",
>> ylab=" body_weight  ( kg )",
>> xaxp=c(0,200,4),
>> yaxp=c(0,100,4),
>> type="p",
>> pch=1,lwd=1.0,
>> cex.lab=1.4, cex.axis=1.2,
>> font.axis=2,
>> cex=1.5,
>> las=1,
>> bty="l",col="firebrick3")
>>
>> boyline<-lm(body_weight ~ body_length, boy)
>> summary(boyline)
>> abline(boyline,col="firebrick3",lwd=2)
>>
>> #~~~
>>   # graph
>> par(mai=c(1.03,1.03,0.4,0.4))
>>
>> par(new=T)
>>
>> plot(girl$body_length, girl$body_weight,
>> xlab=" body_length (cm)",
>> ylab=" body_weight  ( kg )",
>> xaxp=c(0,200,4),
>> yaxp=c(0,100,4),
>> type="p",
>> pch=1,lwd=1.0,
>> cex.lab=1.4, cex.axis=1.2,
>> font.axis=2,
>> cex=1.5,
>> las=1,
>> bty="l",col="saddlebrown")
>>
>>
>> girlline<-lm(body_weight~ body_length, girl)
>> summary(girlline)
>> abline(girlline,col="saddlebrown",lwd=2)
>>
>> __
>> R-help@r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/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] plot xaxp issue

2013-01-07 Thread Elaine Kuo
Thank you Petr, the code is wonderful.

One more question,
you used [as.numeric(together$sex)] to drawing plots many times (Par(new)).
Please kindly advise if there is a similar method to replace drawing
ablines many times.
If not, I am afraid that the ablines will not follow the same Y and
X-axis places.

Elaine


On Mon, Jan 7, 2013 at 8:11 PM, Elaine Kuo  wrote:
> Thanks a lot.
> Please kindly indicate the meaning of the c(8,7).
>
>> together <- rbind(boy, girl)
>> together$sex <- factor(rep(c("boy", "girl"), c(8,7)))
>
> Elaine
>
> On Mon, Jan 7, 2013 at 4:55 PM, PIKAL Petr  wrote:
>> Hi
>>
>> Instead of two plots with par(new = TRUE) try to put boys and girls together 
>> (quite natural thing, they will be pleased 8-)
>>
>> together <- rbind(boy, girl)
>> together$sex <- factor(rep(c("boy", "girl"), c(8,7)))
>>
>> plot(together$body_length, together$body_weight, , 
>> col=c("firebrick3","saddlebrown")[as.numeric(together$sex)], )
>>
>> Regards
>> Petr
>>
>>
>>> -Original Message-
>>> From: r-help-boun...@r-project.org [mailto:r-help-bounces@r-
>>> project.org] On Behalf Of Elaine Kuo
>>> Sent: Monday, January 07, 2013 9:27 AM
>>> To: r-help@r-project.org
>>> Subject: [R] plot xaxp issue
>>>
>>> Hello,
>>>
>>> I have data of Body length and Body weight of 8 boys and 7 girls.
>>>
>>> I want to draw the plot of Body length (for X) and Body weight (for Y)
>>> based on sex.
>>> Then the two plots want to be overlapped for comparison.
>>>
>>> I used the code below but found the unit length of X axis of boy and
>>> girl plot are not the same.
>>> For instance, the length between 0 and 1 of boy plot is larger than
>>> that in girl plot.
>>> The same thing happened to Y axis as well.
>>> (In other words, though axap and yaxp were set to be the same, the
>>> display were not the same.)
>>>
>>> Please kindly advise correction of the code.
>>> Thank you.
>>>
>>> Elaine
>>>
>>> # plot code
>>> boy<-read.csv("H:/boy_data.csv",header=T)
>>> girl<-read.csv("H:/girl_data.csv",header=T)
>>> par(mai=c(1.03,1.03,0.4,0.4))
>>>
>>> plot(boy$body_length, boy$body_weight,
>>> xlab=" body_length (cm)",
>>> ylab=" body_weight  ( kg )",
>>> xaxp=c(0,200,4),
>>> yaxp=c(0,100,4),
>>> type="p",
>>> pch=1,lwd=1.0,
>>> cex.lab=1.4, cex.axis=1.2,
>>> font.axis=2,
>>> cex=1.5,
>>> las=1,
>>> bty="l",col="firebrick3")
>>>
>>> boyline<-lm(body_weight ~ body_length, boy)
>>> summary(boyline)
>>> abline(boyline,col="firebrick3",lwd=2)
>>>
>>> #~~~
>>>   # graph
>>> par(mai=c(1.03,1.03,0.4,0.4))
>>>
>>> par(new=T)
>>>
>>> plot(girl$body_length, girl$body_weight,
>>> xlab=" body_length (cm)",
>>> ylab=" body_weight  ( kg )",
>>> xaxp=c(0,200,4),
>>> yaxp=c(0,100,4),
>>> type="p",
>>> pch=1,lwd=1.0,
>>> cex.lab=1.4, cex.axis=1.2,
>>> font.axis=2,
>>> cex=1.5,
>>> las=1,
>>> bty="l",col="saddlebrown")
>>>
>>>
>>> girlline<-lm(body_weight~ body_length, girl)
>>> summary(girlline)
>>> abline(girlline,col="saddlebrown",lwd=2)
>>>
>>> __
>>> R-help@r-project.org mailing list
>>> https://stat.ethz.ch/mailman/listinfo/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] package metafor: error when setting 'col' and 'at' for a forest plot

2013-01-07 Thread Brian Z Ring
Brilliant, this works very well. Thank you for explaining this so clearly,
all your points are well taken.

Sincerely,
Brian

Brian Z Ring PhD
Professor, Director
Institute of Personalized and Genomic Medicine
College of Life Science
Huazhong University of Science and Technology
Wuhan,  China

-Original Message-
From: Viechtbauer Wolfgang (STAT)
[mailto:wolfgang.viechtba...@maastrichtuniversity.nl] 
Sent: Monday, January 07, 2013 6:16 PM
To: Brian Z Ring; 'David Winsemius'; r-help@r-project.org
Subject: RE: [R] package metafor: error when setting 'col' and 'at' for a
forest plot

Dear Brian,

At the moment, the various forest() functions are not meant to accept 'col'
as an argument. While it is indeed possible to specify a 'col' argument, it
will be passed on via the ... argument to further functions within forest()
and this is where things can go awry. 

To give a more technical explanation, let's look at your example. First,
note that the highest 'at' value you specified is log(10), which is
approximately 2.30. However, the highest upper CI bound in your dataset is
2.9. Therefore, the forest() function wants to draw a right-pointed arrow
for that study using the polygon() function. The code within forest() looks
something like this:

polygon(, col="black", ...)

When you specify a 'col' argument, it is passed on to polygon() via ... --
but col="black" is already hard-coded in forest() and so you get the "formal
argument "col" matched by multiple actual arguments" error.

The issue of changing colors in the forest() functions has come up before:

https://stat.ethz.ch/pipermail/r-help/2011-August/287788.html

In fact, the problems associated with allowing the user to specify colors in
high-level plotting functions, such as forest(), was nicely discussed in a
recent article by Paul Murrell in the R Journal:

http://journal.r-project.org/archive/2012-2/RJournal_2012-2_Murrell.pdf

I may add functionality to the forest() functions to handle a 'col' argument
for specifying row colors, but I will have to evaulate carefully how well
that will work. Also, I would then need to change the current 'col' argument
in the forest.rma() function -- breaking backwards compatability :/

Back to your example. You can actually get it to work if you use:

forest(forest$OR, ci.lb=forest$Low, ci.ub=forest$High,  at=log(c(.05, .25,
1, 20)), slab=forest$SNP, atransf=exp)

in which case the values of 'at' encompass the CI bounds of all studies.
Then the forest() function does not use polygon(), but segments() and here
'col' is not hard-coded. One (probably) unintended consequence of using
'col' this way is that the x axis label is also given a different color. 

A few other notes:

1) Calling your data frame 'forest' is probably not a good idea -- since
that is the name of a function.

2) Your data appear to be (raw -- not log transformed) odds ratios and the
corresponding CI bounds. So, you should not be using the atransf arguments
-- or you would be transforming your odds ratios to exp(ORs). Probably
closer to what you want is:

dat <- read.table(header=TRUE, text="
SNP   Group   High   Low   OR
rs1137101   A   1.21   0.87   1.03
rs1137101   B   2.11   1.21   1.6
rs1137101   C   2.9   1.42   2.03
rs1042522   A   1.12   0.84   0.97
rs1042522   B   1.15   0.79   0.95
rs1042522   C   0.92   0.5   0.7
rs1625895   A   1.14   0.76   0.93
rs1625895   B   1.15   0.75   0.93
rs1625895   C   NA   NA   NA
ACEI/D   A   1.55   0.79   1.11
ACEI/D   B   1.25   0.76   0.98
ACEI/D   C   0.85   0.41   0.59
")

forest(dat$OR, ci.lb=dat$Low, ci.ub=dat$High, col=c(1,2,3), at=c(.25, 1, 2,
3), slab=dat$SNP, refline=1, xlab=" ", xlim=c(-1.5,5.5)) mtext("Odds Ratio",
side=1, line=2.5, adj=0.45)

where I manually add the x axis label, so that it is in black. Again,
though, using 'col' this way is technically not intended (even though it
does work here).

I hope this helps!

Best,
Wolfgang

--   
Wolfgang Viechtbauer, Ph.D., Statistician   
Department of Psychiatry and Psychology   
School for Mental Health and Neuroscience   
Faculty of Health, Medicine, and Life Sciences   
Maastricht University, P.O. Box 616 (VIJV1)   
6200 MD Maastricht, The Netherlands   
+31 (43) 388-4170 | http://www.wvbauer.com   

> -Original Message-
> From: r-help-boun...@r-project.org 
> [mailto:r-help-boun...@r-project.org]
> On Behalf Of Brian Z Ring
> Sent: Sunday, January 06, 2013 05:48
> To: 'David Winsemius'; r-help@r-project.org
> Subject: Re: [R] package metafor: error when setting 'col' and 'at' 
> for a forest plot
> 
> Good suggestion, and yes, the given examples work. Using a 
> modification of one of the supplied examples, where in the final line 
> I only added the parameter col=c(1,2,3):
> 
> data(dat.bcg)
> dat <- escalc(measure="RR", ai=tpos, bi=tneg, ci=cpos, di=cneg, 
> data=dat.bcg, append=TRUE) forest(dat$yi, dat$vi, 
> slab=paste(dat$author, dat$year, sep=", "), atransf=exp, 
> at=log(c(.05,.25,1,4,20)), col=c(1,2,3))
> 
> I get a plot sim

Re: [R] plot xaxp issue

2013-01-07 Thread PIKAL Petr
Hi

> -Original Message-
> From: Elaine Kuo [mailto:elaine.kuo...@gmail.com]
> Sent: Monday, January 07, 2013 1:45 PM
> To: PIKAL Petr
> Cc: r-help@r-project.org
> Subject: Re: [R] plot xaxp issue
> 
> Thank you Petr, the code is wonderful.
> 
> One more question,
> you used [as.numeric(together$sex)] to drawing plots many times
> (Par(new)).

It is not about drawing plot many times but coding points or graphic objects by 
some factor. In your case sex.

> Please kindly advise if there is a similar method to replace drawing
> ablines many times.

Instead of
boyline<-lm(body_weight ~ body_length, boy)

use collective data frame together and subset only one sex.

boyline<-lm(body_weight ~ body_length, data=together, subset(sex=="boy"))
abline(boyline,col="firebrick3",lwd=2)

In case of more than two sexes (some SF authors mentioned it is possible) you 
can use simple loop.

for (i in 1:2)) {
subs <- together$sex==levels(together$sex)[i]
line<-lm(body_weight ~ body_length, data=together, subset=subs)
abline(line, col=c("firebrick3","saddlebrown")[i],lwd=2) 
}

> If not, I am afraid that the ablines will not follow the same Y and X-
> axis places.

They will. You can try it yourself. Sometimes R functions are quite close to 
magic. 

Regards
Petr

> 
> Elaine
> 
> 
> On Mon, Jan 7, 2013 at 8:11 PM, Elaine Kuo 
> wrote:
> > Thanks a lot.
> > Please kindly indicate the meaning of the c(8,7).
> >
> >> together <- rbind(boy, girl)
> >> together$sex <- factor(rep(c("boy", "girl"), c(8,7)))
> >
> > Elaine
> >
> > On Mon, Jan 7, 2013 at 4:55 PM, PIKAL Petr 
> wrote:
> >> Hi
> >>
> >> Instead of two plots with par(new = TRUE) try to put boys and girls
> >> together (quite natural thing, they will be pleased 8-)
> >>
> >> together <- rbind(boy, girl)
> >> together$sex <- factor(rep(c("boy", "girl"), c(8,7)))
> >>
> >> plot(together$body_length, together$body_weight, ,
> >> col=c("firebrick3","saddlebrown")[as.numeric(together$sex)], )
> >>
> >> Regards
> >> Petr
> >>
> >>
> >>> -Original Message-
> >>> From: r-help-boun...@r-project.org [mailto:r-help-bounces@r-
> >>> project.org] On Behalf Of Elaine Kuo
> >>> Sent: Monday, January 07, 2013 9:27 AM
> >>> To: r-help@r-project.org
> >>> Subject: [R] plot xaxp issue
> >>>
> >>> Hello,
> >>>
> >>> I have data of Body length and Body weight of 8 boys and 7 girls.
> >>>
> >>> I want to draw the plot of Body length (for X) and Body weight (for
> >>> Y) based on sex.
> >>> Then the two plots want to be overlapped for comparison.
> >>>
> >>> I used the code below but found the unit length of X axis of boy
> and
> >>> girl plot are not the same.
> >>> For instance, the length between 0 and 1 of boy plot is larger than
> >>> that in girl plot.
> >>> The same thing happened to Y axis as well.
> >>> (In other words, though axap and yaxp were set to be the same, the
> >>> display were not the same.)
> >>>
> >>> Please kindly advise correction of the code.
> >>> Thank you.
> >>>
> >>> Elaine
> >>>
> >>> # plot code
> >>> boy<-read.csv("H:/boy_data.csv",header=T)
> >>> girl<-read.csv("H:/girl_data.csv",header=T)
> >>> par(mai=c(1.03,1.03,0.4,0.4))
> >>>
> >>> plot(boy$body_length, boy$body_weight,
> >>> xlab=" body_length (cm)",
> >>> ylab=" body_weight  ( kg )",
> >>> xaxp=c(0,200,4),
> >>> yaxp=c(0,100,4),
> >>> type="p",
> >>> pch=1,lwd=1.0,
> >>> cex.lab=1.4, cex.axis=1.2,
> >>> font.axis=2,
> >>> cex=1.5,
> >>> las=1,
> >>> bty="l",col="firebrick3")
> >>>
> >>> boyline<-lm(body_weight ~ body_length, boy)
> >>> summary(boyline)
> >>> abline(boyline,col="firebrick3",lwd=2)
> >>>
> >>>
> #~~~
> >>>   # graph
> >>> par(mai=c(1.03,1.03,0.4,0.4))
> >>>
> >>> par(new=T)
> >>>
> >>> plot(girl$body_length, girl$body_weight,
> >>> xlab=" body_length (cm)",
> >>> ylab=" body_weight  ( kg )",
> >>> xaxp=c(0,200,4),
> >>> yaxp=c(0,100,4),
> >>> type="p",
> >>> pch=1,lwd=1.0,
> >>> cex.lab=1.4, cex.axis=1.2,
> >>> font.axis=2,
> >>> cex=1.5,
> >>> las=1,
> >>> bty="l",col="saddlebrown")
> >>>
> >>>
> >>> girlline<-lm(body_weight~ body_length, girl)
> >>> summary(girlline)
> >>> abline(girlline,col="saddlebrown",lwd=2)
> >>>
> >>> __
> >>> R-help@r-project.org mailing list
> >>> https://stat.ethz.ch/mailman/listinfo/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] Function to tabulate t-test (1) How can I get dependent variable name (2) How can I make data frame know to function

2013-01-07 Thread John Sorkin
Windows 7
R 2.12.1
I am trying to write a function (see sample code below) that will take the 
output of a t-test and produce results suitable for a table.
I have two questions
(1) You will note that the name of the outcome variable, which is "value" in 
the input is replaced by the string "outcome by class" in the data frame 
produced by my function. How can I make my function put the name of the 
variable being analyzed,, i.e "value" in the output data frame
(2) How can I pass the entire input data to the function so my call to the 
function will not have to be in its current ugly form, i.e.  
Table1(data$value,data$sex)
and instead could just be
Table1(value,sex,data=data)
 
 
x <- data.frame(value=rnorm(20)  ,sex=rep("Male",  20))
y <- data.frame(value=rnorm(20,4),sex=rep("Female",20))
data <- rbind(x,y)
temp <- t.test(value~sex,data=data)
temp
v<-data.frame(dep=temp$data.name,
   female=temp$estimate[1],male=temp$estimate[2],
   p=temp$p.value,
   CILow=temp$conf.int[1],CIHigh=temp$conf.int[2])
row.names(v) <- NULL
 
 
 
Table1 <- function(outcome,class) {
temp <- t.test(outcome~ class)
mydf <- data.frame(dep=temp$data.name,
 female=temp$estimate[1],male=temp$estimate[2],
 p=temp$p.value,
 CILow=temp$conf.int[1],CIHigh=temp$conf.int[2])
row.names(mydf) <- NULL
mydf}
Table1(data$value,data$sex)
 
 
 
Thank you,
John
 
 
 
John David Sorkin M.D., Ph.D.
Chief, Biostatistics and Informatics
University of Maryland School of Medicine Division of Gerontology
Baltimore VA Medical Center
10 North Greene Street
GRECC (BT/18/GR)
Baltimore, MD 21201-1524
(Phone) 410-605-7119
(Fax) 410-605-7913 (Please call phone number above prior to faxing)
Confidentiality Statement:
This email message, including any attachments, is for the sole use of the 
intended recipient(s) and may contain confidential and privileged information.  
Any unauthorized use, disclosure or distribution is prohibited.  If you are not 
the intended recipient, please contact the sender by reply email and destroy 
all copies of the original message. 
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] checker/chequer board pattern as a colour

2013-01-07 Thread Federico Calboli
Hi All,

is there a reasonably simple way of using a black and white chequer/checker 
board pattern as a colour:

barplot(mydata, col = c('red', 'blue' 'checkerboard'))

?

BW

F

--
Federico C. F. Calboli
Neuroepidemiology and Ageing Research
Imperial College, St. Mary's Campus
Norfolk Place, London W2 1PG

Tel +44 (0)20 75941602   Fax +44 (0)20 75943193

f.calboli [.a.t] imperial.ac.uk
f.calboli [.a.t] 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] Function to tabulate t-test (1) How can I get dependent variable name (2) How can I make data frame know to function

2013-01-07 Thread Bert Gunter
John:

R2.12 ??? Time to update.

1. ?t.test and note the last entry in the Value section.

2.  ?with

-- Bert

On Mon, Jan 7, 2013 at 6:44 AM, John Sorkin  wrote:
> Windows 7
> R 2.12.1
> I am trying to write a function (see sample code below) that will take the 
> output of a t-test and produce results suitable for a table.
> I have two questions
> (1) You will note that the name of the outcome variable, which is "value" in 
> the input is replaced by the string "outcome by class" in the data frame 
> produced by my function. How can I make my function put the name of the 
> variable being analyzed,, i.e "value" in the output data frame
> (2) How can I pass the entire input data to the function so my call to the 
> function will not have to be in its current ugly form, i.e.
> Table1(data$value,data$sex)
> and instead could just be
> Table1(value,sex,data=data)
>
>
> x <- data.frame(value=rnorm(20)  ,sex=rep("Male",  20))
> y <- data.frame(value=rnorm(20,4),sex=rep("Female",20))
> data <- rbind(x,y)
> temp <- t.test(value~sex,data=data)
> temp
> v<-data.frame(dep=temp$data.name,
>female=temp$estimate[1],male=temp$estimate[2],
>p=temp$p.value,
>CILow=temp$conf.int[1],CIHigh=temp$conf.int[2])
> row.names(v) <- NULL
>
>
>
> Table1 <- function(outcome,class) {
> temp <- t.test(outcome~ class)
> mydf <- data.frame(dep=temp$data.name,
>  female=temp$estimate[1],male=temp$estimate[2],
>  p=temp$p.value,
>  CILow=temp$conf.int[1],CIHigh=temp$conf.int[2])
> row.names(mydf) <- NULL
> mydf}
> Table1(data$value,data$sex)
>
>
>
> Thank you,
> John
>
>
>
> John David Sorkin M.D., Ph.D.
> Chief, Biostatistics and Informatics
> University of Maryland School of Medicine Division of Gerontology
> Baltimore VA Medical Center
> 10 North Greene Street
> GRECC (BT/18/GR)
> Baltimore, MD 21201-1524
> (Phone) 410-605-7119
> (Fax) 410-605-7913 (Please call phone number above prior to faxing)
> Confidentiality Statement:
> This email message, including any attachments, is for ...{{dropped:27}}

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Function to tabulate t-test (1) How can I get dependent variable name (2) How can I make data frame know to function

2013-01-07 Thread Rui Barradas
Hello,

A simple change to your function partly answers to both questions.


Table2 <- function(formula, data) {
 temp <- t.test(formula, data)
 mydf <- data.frame(dep=temp$data.name,
  female=temp$estimate[1],male=temp$estimate[2],
  p=temp$p.value,
  CILow=temp$conf.int[1],CIHigh=temp$conf.int[2])
 row.names(mydf) <- NULL
 mydf
}
Table2(value ~ sex, data)


If you want the outcome variable to be called just "value", inside the 
function use

dname <- rownames(attr(terms(formula), "factors"))[1]
mydf <- data.frame(dep=dname, [...etc...]


Hope this helps,

Rui Barradas


Em 07-01-2013 14:44, John Sorkin escreveu:
> Windows 7
> R 2.12.1
> I am trying to write a function (see sample code below) that will take the 
> output of a t-test and produce results suitable for a table.
> I have two questions
> (1) You will note that the name of the outcome variable, which is "value" in 
> the input is replaced by the string "outcome by class" in the data frame 
> produced by my function. How can I make my function put the name of the 
> variable being analyzed,, i.e "value" in the output data frame
> (2) How can I pass the entire input data to the function so my call to the 
> function will not have to be in its current ugly form, i.e.
> Table1(data$value,data$sex)
> and instead could just be
> Table1(value,sex,data=data)
>   
>   
> x <- data.frame(value=rnorm(20)  ,sex=rep("Male",  20))
> y <- data.frame(value=rnorm(20,4),sex=rep("Female",20))
> data <- rbind(x,y)
> temp <- t.test(value~sex,data=data)
> temp
> v<-data.frame(dep=temp$data.name,
> female=temp$estimate[1],male=temp$estimate[2],
> p=temp$p.value,
> CILow=temp$conf.int[1],CIHigh=temp$conf.int[2])
> row.names(v) <- NULL
>   
>   
>   
> Table1 <- function(outcome,class) {
> temp <- t.test(outcome~ class)
> mydf <- data.frame(dep=temp$data.name,
>   female=temp$estimate[1],male=temp$estimate[2],
>   p=temp$p.value,
>   CILow=temp$conf.int[1],CIHigh=temp$conf.int[2])
> row.names(mydf) <- NULL
> mydf}
> Table1(data$value,data$sex)
>   
>   
>   
> Thank you,
> John
>   
>   
>   
> John David Sorkin M.D., Ph.D.
> Chief, Biostatistics and Informatics
> University of Maryland School of Medicine Division of Gerontology
> Baltimore VA Medical Center
> 10 North Greene Street
> GRECC (BT/18/GR)
> Baltimore, MD 21201-1524
> (Phone) 410-605-7119
> (Fax) 410-605-7913 (Please call phone number above prior to faxing)
> Confidentiality Statement:
> This email message, including any attachments, is for ...{{dropped:13}}

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Function to tabulate t-test (1) How can I get dependent variable name (2) How can I make data frame know to function

2013-01-07 Thread Rui Barradas
I've just realized that you're swapping female and male in the creation 
of the results data frame.

It should be

Table2 <- function(formula, data) {
dname <- rownames(attr(terms(formula), "factors"))[1]
temp <- t.test(formula, data)
mydf <- data.frame(dep=temp$data.name,
 female=temp$estimate[2], male=temp$estimate[1],
 p=temp$p.value,
 CILow=temp$conf.int[1],CIHigh=temp$conf.int[2])
row.names(mydf) <- NULL
mydf
}


Rui Barradas

Em 07-01-2013 15:39, Rui Barradas escreveu:

Hello,

A simple change to your function partly answers to both questions.


Table2 <- function(formula, data) {
  temp <- t.test(formula, data)
  mydf <- data.frame(dep=temp$data.name,
   female=temp$estimate[1],male=temp$estimate[2],
   p=temp$p.value,
   CILow=temp$conf.int[1],CIHigh=temp$conf.int[2])
  row.names(mydf) <- NULL
  mydf
}
Table2(value ~ sex, data)


If you want the outcome variable to be called just "value", inside the
function use

dname <- rownames(attr(terms(formula), "factors"))[1]
mydf <- data.frame(dep=dname, [...etc...]


Hope this helps,

Rui Barradas


Em 07-01-2013 14:44, John Sorkin escreveu:

Windows 7
R 2.12.1
I am trying to write a function (see sample code below) that will take the 
output of a t-test and produce results suitable for a table.
I have two questions
(1) You will note that the name of the outcome variable, which is "value" in the input is replaced 
by the string "outcome by class" in the data frame produced by my function. How can I make my 
function put the name of the variable being analyzed,, i.e "value" in the output data frame
(2) How can I pass the entire input data to the function so my call to the 
function will not have to be in its current ugly form, i.e.
Table1(data$value,data$sex)
and instead could just be
Table1(value,sex,data=data)
   
   
x <- data.frame(value=rnorm(20)  ,sex=rep("Male",  20))

y <- data.frame(value=rnorm(20,4),sex=rep("Female",20))
data <- rbind(x,y)
temp <- t.test(value~sex,data=data)
temp
v<-data.frame(dep=temp$data.name,
 female=temp$estimate[1],male=temp$estimate[2],
 p=temp$p.value,
 CILow=temp$conf.int[1],CIHigh=temp$conf.int[2])
row.names(v) <- NULL
   
   
   
Table1 <- function(outcome,class) {

temp <- t.test(outcome~ class)
mydf <- data.frame(dep=temp$data.name,
   female=temp$estimate[1],male=temp$estimate[2],
   p=temp$p.value,
   CILow=temp$conf.int[1],CIHigh=temp$conf.int[2])
row.names(mydf) <- NULL
mydf}
Table1(data$value,data$sex)
   
   
   
Thank you,

John
   
   
   
John David Sorkin M.D., Ph.D.

Chief, Biostatistics and Informatics
University of Maryland School of Medicine Division of Gerontology
Baltimore VA Medical Center
10 North Greene Street
GRECC (BT/18/GR)
Baltimore, MD 21201-1524
(Phone) 410-605-7119
(Fax) 410-605-7913 (Please call phone number above prior to faxing)
Confidentiality Statement:
This email message, including any attachments, is for ...{{dropped:13}}

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


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


Re: [R] how to aggregate T-test result in an elegant way?

2013-01-07 Thread Yao He
Hi,arun

Yes , I just want to do the t.test
I think maybe  it is not necessary to generate a 3D array from the raw
data.frame by acast() at first

Thanks a lot

2013/1/7 arun :
> Hi Yao,
>
> It's okay.
>
> How did you generate the 3 D array?
> Using ?acast()
>
> I am not sure I understand your question "
>
> if you meet a t-test task as I described  , is that generate a
> high-dimension array  a good way ?"
>
> Do you want to do the t-test in the melt dataset?
>
> b<- read.table(text="
> IDO2variablevalue
> 1TWF2H513% EW.INCU49.38
> 2TWF2H613% EW.INCU48.02
> 3TWF2H1913%EW.INCU51.44
> 280TWF2H10113% EW.17.542.26
> 281TWF2H10513%EW.17.543.52
> 282TWF2H10613% EW.17.542.83
> 472TWF2N10221% EW.17.545.97
> 473TWF2N10421%EW.17.5 43.32
> 474TWF2N10621% EW.17.548.63
> 689TWF2N221% EMW19.57
> 690TWF2N621%EMW18.07
> 691TWF2N1021%EMW15.4
> 491TWF2H513%EMW15.61
> 492TWF2H613%EMW13.41
> 493TWF2H1913%EMW14.03
> 199TWF2N221%EW.INCU48.69
> 200TWF2N621%EW.INCU50.52
> 201TWF2N1021%EW.INCU42.04
> ",sep="",header=TRUE,stringsAsFactors=FALSE)
>  res<-lapply(lapply(split(b,b$variable),function(x) 
> t.test(x$value[x$O2=="13%"],x$value[x$O2=="21%"])),function(x) 
> data.frame(mean=x$estimate,p.value=x$p.value))
> res1<-do.call(rbind,res)
> row.names(res1)[grep("mean of 
> x",row.names(res1))]<-gsub("(.*\\.).*$","\\113%",row.names(res1)[grep("mean 
> of x",row.names(res1))])
> row.names(res1)[grep("mean of 
> y",row.names(res1))]<-gsub("(.*\\.).*$","\\121%",row.names(res1)[grep("mean 
> of y",row.names(res1))])
> res1
> #meanp.value
> #EMW.13% 14.35000 0.09355374
> #EMW.21% 17.68000 0.09355374
> #EW.17.5.13% 42.87000 0.17464018
> #EW.17.5.21% 45.97333 0.17464018
> #EW.INCU.13% 49.61333 0.43689727
> #EW.INCU.21% 47.08333 0.43689727
>
> A.K.
>
>
>
> - Original Message -
> From: Yao He 
> To: arun 
> Cc: R help 
> Sent: Monday, January 7, 2013 4:00 AM
> Subject: Re: [R] how to aggregate T-test result in an elegant way?
>
> Hi, arun
> I'm so sorry for that isn't helpful.
> One of question is that I don't know how  to subset a small part as it
> is a 3-dimension array so I just show the structure of that.
> I tried  dput()  to a file , then what should I do for subsetting it?
>
> Another question is :
> My rawdata is a "melt" dataframe like that:
> IIDO2variablevalue
> 1TWF2H513% EW.INCU49.38
> 2TWF2H613% EW.INCU48.02
> 3TWF2H1913% EW.INCU51.44
> 280TWF2H10113% EW.17.542.26
> 281TWF2H10513% EW.17.5 43.52
> 282TWF2H10613% EW.17.542.83
> 472TWF2N10221% EW.17.545.97
> 473TWF2N10421% EW.17.5 43.32
> 474TWF2N10621% EW.17.548.63
> 689TWF2N221%  EMW19.57
> 690TWF2N621% EMW18.07
> 691TWF2N1021%EMW15.4
> 491TWF2H5 13%EMW15.61
> 492TWF2H613%EMW13.41
> 493TWF2H1913%EMW14.03
> 199TWF2N221%EW.INCU48.69
> 200TWF2N621%EW.INCU50.52
> 201TWF2N1021%EW.INCU42.04
>
> if you meet a t-test task as I described  , is that generate a
> high-dimension array  a good way ?
> Thank you!
>
> Yao He
> 2013/1/7 arun :
>> HI,
>> I tried to create an example dataset (as you didn't provide the data).
>> set.seed(25)
>> a<-array(sample(1:50,60,replace=TRUE),dim=c(2,10,3))
>> dimnames(a)[[1]]<-c("13%","21%")
>> dimnames(a)[[2]]<-paste("TWF2H",101:110,sep="")
>> dimnames(a)[[3]]<-c("EW.INCU","EW.17.5","EMW")
>>
>>
>> str(a)
>> # int [1:2, 1:10, 1:3] 21 35 8 45 7 50 32 17 4 15 ...
>>  #- attr(*, "dimnames")=List of 3
>>   #..$ : chr [1:2] "13%" "21%"
>>   .#.$ : chr [1:10] "TWF2H101" "TWF2H102" "TWF2H103" "TWF2H104" ...
>>   #..$ : chr [1:3] "EW.INCU" "EW.17.5" "EMW"
>>
>> res<-lapply(lapply(seq_len(dim(a)[3]),function(i) 
>> t.test(a[dimnames(a)[[1]][1],,i],a[dimnames(a)[[1]][2],,i])),function(x) 
>> data.frame(mean=x$estimate,p.value=x$p.value))
>> res1<-do.call(rbind,res)
>>   row.names(res1)[grep("mean of 
>> x",row.names(res1))]<-gsub("(.*\\.).*$","\\113%",row.names(res1)[grep("mean 
>> of x",row.names(res1))])
>>  row.names(res1)[grep("mean of 
>> y",row.names(res1))]<-gsub("(.*\\.).*$","\\121%",row.names(res1)[grep("mean 
>> of y",row.names(res1))])
>> res1
>> #mean   p.value
>> #EW.INCU.13% 22.3 0.2754842
>> #EW.INCU.21% 29.3 0.2754842
>> #EW.17.5.13% 20.5 0.4705772
>> #EW.17.5.21% 16.0 0.4705772
>> #EMW.13% 23.9 0.9638679
>> #EMW.21% 24.2 0.9638679
>> A.K.
>>
>>
>>
>>
>> - Original Message -
>> From: Yao He 
>> To: arun 
>> Cc: R help 
>> Sent: Sunday, January 6, 2013 11:21 PM
>> Subject: Re: [R] h

[R] list of lists to matrix

2013-01-07 Thread eliza botto

dear R family,
[a text file has been attached for better understanding]
i have a list of 16 and each of of that is further subdivided into variable 
number of lists. So, i have a kind of list into lists phenomenon.
[[1]]$'1'
1   2  3   4   5  6
7  8   9

[[1]]$'2'
1   2  3   4   5  6
7  8   9
i want to convert both these sublists into one column and then cbind it in the 
following way
col1 col2
1   1
2   2
3   3..
9 9

i want to the same operations on all the 16 lists. 
thanks in advance,

elisa dear R family,

i have a list of 16 and each of of that is further subdivided into variable 
number of lists. So, i have a kind of list of lists phenomenon.

[[1]]$'1'

1   2  3
   
4   5  6

7  8   9


[[1]]$'2'

1   2  3
   
4   5  6

7  8   9

i want to convert both these sublists into one column and then cbind it in the 
following way

col1 col2

1   1

2   2

3   3
.
.

9 9


i want to the same operations on all the 16 lists. 

thanks in advance,

elisa__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Have problem to do loop to generate transformed chi-squared variates

2013-01-07 Thread Agnes Ayang
Hello R-helpers,

I need to generate standard variates normal to 'create' chi-squared variates. 
To make you more understand,

(1)   a<-rnorm(3,0,1)

*after do (1), I need to squared and summed the three values. My problem is, 
how am I going to continue the programming if I had to repeat the process for 
15 times, which in the end I will get 15 values from the whole programme.Hope 
you can help. 


Sincerely,
Agnes
[[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] Have problem to do loop to generate transformed chi-squared variates

2013-01-07 Thread R. Michael Weylandt
Hi Agnes,

I think this is likely homework, right? If that's the case, we're
really not supposed to give you help.

If not, why not simply use rchisq to generate chi-sq variates exactly?

MW

On Mon, Jan 7, 2013 at 4:13 PM, Agnes Ayang  wrote:
> Hello R-helpers,
>
> I need to generate standard variates normal to 'create' chi-squared variates. 
> To make you more understand,
>
> (1)   a<-rnorm(3,0,1)
>
> *after do (1), I need to squared and summed the three values. My problem is, 
> how am I going to continue the programming if I had to repeat the process for 
> 15 times, which in the end I will get 15 values from the whole programme.Hope 
> you can help.
>
>
> Sincerely,
> Agnes
> [[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] renumber a list of numbers

2013-01-07 Thread Charles Determan Jr
Greetings R users,

I am trying to renumber my groups within the file shown below.  The groups
are currently set as 8,9,10,etc.  I would like to renumber this as
1,2,3,etc.  I have searched the help files and only come across using the
rownames to renumber the values but I need to match values.  Any assistance
is always appreciated,

Regards,
Charles

structure(list(Group = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L,
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L,
1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L), .Label = c("NO", "YES"
), class = "factor"), Event_name = c(8, 9, 10, 11, 12, 13, 14,
15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30,
31, 33, 34, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20,
21, 22, 23, 24, 25, 26, 27, 28, 29, 30), Glucose.n = c(26, 3,
0, 26, 1, 25, 26, 25, 25, 26, 26, 25, 23, 23, 24, 24, 26, 26,
25, 25, 26, 26, 24, 21, 7, 12, 4, 0, 0, 4, 0, 4, 4, 4, 4, 4,
4, 4, 3, 4, 4, 4, 3, 3, 2, 2, 2, 1, 1), Glucose.m = c(92.5,
90.3,
NaN, 97.2307692307692, 116, 97.84, 107.653846153846, 105.32,
102.6, 94.6538461538462, 96.076923076923, 92.24, 87.5652173913043,
79.3913043478261, 81.29167, 77.5, 75.9230769230769,
74.4615384615385,
72.68, 76.32, 74.9615384615385, 72.2307692307692, 92.54167,
105.619047619048, 93.4285714285714, 96.5, 90, NaN, NaN, 86.5,
NaN, 87, 90.25, 92.5, 98.75, 95.75, 94, 88.25, 54.3,
52, 74.5, 77.75, 81, 97.3, 82.5, 85, 66.5, 51, 81
), Glucose.sd = c(18.9256439784753, 27.5922694487665, NA, 25.3050314242961,
NA, 17.3917605012642, 21.027491163127, 12.0094407308029, 28.0728219695373,
17.7334538264655, 10.7700439253443, 12.7778454104490, 11.0075432274935,
14.6992242542214, 12.2739709270814, 10.9266328000819, 10.4457573279225,
13.1338669682033, 8.2194890352138, 19.9556174213344, 17.6079090620795,
10.9299869800753, 19.3052801217623, 29.7883806046522, 17.2032665779607,
18.4563366797521, 18.0554700852678, NA, NA, 19.7399763593239,
NA, 22.3159136044214, 28.5116935075184, 21.7638844572072, 12.2848144742469,
12.0933866224478, 15.5777619273972, 11.842719282327, 39.1066916694999,
32.0936130717624, 66.8755062286136, 53.7796429887741, 17.6918060129541,
11.5902257671425, 34.6482322781408, 9.89949493661167, 26.1629509039023,
NA, NA), Glucose.se = c(3.71162415205983, 15.9304041937980, NA,
4.96272496248326, NA, 3.47835210025284, 4.12383029856479, 2.40188814616057,
5.61456439390745, 3.47781642709786, 2.11217938990701, 2.6908208981,
2.29523142628873, 3.06500013246251, 2.50541382409208, 2.23038958057991,
2.04858155574352, 2.57576322922309, 1.64389780704276, 3.99112348426689,
3.45319507311964, 2.14354680364304, 3.94067380331773, 6.50035756909396,
6.50222358617618, 5.32788547515463, 9.0277350426339, NA, NA,
9.86998817966195, NA, 11.1579568022107, 14.2558467537592, 10.8819422286036,
6.14240723712346, 6.04669331122391, 7.7096369861, 5.92135964116351,
22.5782589625015, 16.0468065358812, 33.4377531143068, 26.8898214943871,
10.2143689640297, 6.69161996662824, 24.5, 7, 18.5, NA, NA)), .Names =
c("Group",
"Event_name", "Glucose.n", "Glucose.m", "Glucose.sd", "Glucose.se"
), row.names = c(9L, 10L, 11L, 12L, 13L, 14L, 15L, 16L, 17L,
18L, 19L, 20L, 21L, 22L, 23L, 24L, 25L, 26L, 27L, 28L, 29L, 30L,
31L, 32L, 33L, 34L, 41L, 42L, 43L, 44L, 45L, 46L, 47L, 48L, 49L,
50L, 51L, 52L, 53L, 54L, 55L, 56L, 57L, 58L, 59L, 60L, 61L, 62L,
63L), class = "data.frame")

[[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] list of lists to matrix

2013-01-07 Thread R. Michael Weylandt
Something like

do.call(cbind, lists)

?

MW

On Mon, Jan 7, 2013 at 4:13 PM, eliza botto  wrote:
>
> dear R family,
> [a text file has been attached for better understanding]
> i have a list of 16 and each of of that is further subdivided into variable 
> number of lists. So, i have a kind of list into lists phenomenon.
> [[1]]$'1'
> 1   2  3   4   5  6
> 7  8   9
>
> [[1]]$'2'
> 1   2  3   4   5  6
> 7  8   9
> i want to convert both these sublists into one column and then cbind it in 
> the following way
> col1 col2
> 1   1
> 2   2
> 3   3..
> 9 9
>
> i want to the same operations on all the 16 lists.
> thanks in advance,
>
> elisa
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/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] list of lists to matrix

2013-01-07 Thread eliza botto

Thanks arun and weylandt,it perfectly worked out..
elisa

> From: michael.weyla...@gmail.com
> Date: Mon, 7 Jan 2013 16:37:53 +
> Subject: Re: [R] list of lists to matrix
> To: eliza_bo...@hotmail.com
> CC: r-help@r-project.org
> 
> Something like
> 
> do.call(cbind, lists)
> 
> ?
> 
> MW
> 
> On Mon, Jan 7, 2013 at 4:13 PM, eliza botto  wrote:
> >
> > dear R family,
> > [a text file has been attached for better understanding]
> > i have a list of 16 and each of of that is further subdivided into variable 
> > number of lists. So, i have a kind of list into lists phenomenon.
> > [[1]]$'1'
> > 1   2  3   4   5  6
> > 7  8   9
> >
> > [[1]]$'2'
> > 1   2  3   4   5  6
> > 7  8   9
> > i want to convert both these sublists into one column and then cbind it in 
> > the following way
> > col1 col2
> > 1   1
> > 2   2
> > 3   3..
> > 9 9
> >
> > i want to the same operations on all the 16 lists.
> > thanks in advance,
> >
> > elisa
> > __
> > R-help@r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/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] random effects model

2013-01-07 Thread rex2013
Hi A.K

Below is the comment I get, not sure why.

BP.sub3 is the stacked data without the missing values.

BP.geese3 <- geese(HiBP~time*MaternalAge,data=BP.sub3,id=CODEA,
family=binomial, corstr="unstructured", na.action=na.omit)Error in
`contrasts<-`(`*tmp*`, value = contr.funs[1 + isOF[nn]]) :
  contrasts can be applied only to factors with 2 or more levels

Even though age has 3 levels; time has 14 years & 21 years; HIBP is a
binary response outcome.

2) When you mentioned summary(m1)$mean["p"] what did the p mean? i
used this in one of the gee command, it produced NA as answer?

Many thanks



On Mon, Jan 7, 2013 at 5:26 AM, arun kirshna [via R] <
ml-node+s789695n4654795...@n4.nabble.com> wrote:

> Hi,
>
> I am  not very familiar with the geese/geeglm().  Is it from
> library(geepack)?
> Regarding your question:
> "
> Can you tell me if I can use the geese or geeglm function with this data
> eg: : HIBP~ time* Age
> Here age is a factor with 3 levels, time: 2 levels, HIBP = yes/no.
>
> From your original data:
> BP_2b<-read.csv("BP_2b.csv",sep="\t")
> head(BP_2b,2)
> #  CODEA Sex MaternalAge Education Birthplace AggScore IntScore Obese14
> #1 1  NA   3 4  1   NA   NA  NA
> #2 3   2   3 3  100   0
>  # Overweight14 Overweight21 Obese21 hibp14 hibp21
> #1   NA   NA  NA NA NA
> #201   0  0  0
>
> If I understand your new classification:
> BP.stacknormal<- subset(BP_2b,Obese14==0 & Overweight14==0 & Obese21==0 &
> Overweight21==0)
> BP.stackObese <- subset(BP_2b,(Obese14==1& Overweight14==0 &
> Obese14==1&Overweight14==1)|(Obese14==1&Overweight14==1 & Obese21==1 &
> Overweight21==0)|(Obese14==1&Overweight14==0 & Obese21==0 &
> Overweight21==0)|(Obese14==0 & Overweight14==0 & Obese21==1 &
> Overweight21==0)|(Obese14==0 & Overweight14==0 & Obese21==1 &
> Overweight21==1)|(Obese14==0 & Overweight14==1 & Obese21==1
> &Overweight21==1)|(Obese14==1& Overweight14==1 & Obese21==1&
> Overweight21==1)) #check whether there are more classification that fits to
> #Obese
>  BP.stackOverweight <- subset(BP_2b,(Obese14==0 & Overweight14==1 &
> Obese21==0 & Overweight21==1)|(Obese14==0 &Overweight14==1 & Obese21==0 &
> Overweight21==0)|(Obese14==0 & Overweight14==0 & Obese21==0 &
> Overweight21==1))
> BP.stacknormal$Categ<-"Normal"
> BP.stackObese$Categ<-"Obese"
> BP.stackOverweight$Categ <- "Overweight"
>  
> BP.newObeseOverweightNormal<-na.omit(rbind(BP.stacknormal,BP.stackObese,BP.stackOverweight))
>
>  nrow(BP.newObeseOverweightNormal)
> #[1] 1581
> BP.stack3 <-
> reshape(BP.newObeseOverweightNormal,idvar="CODEA",timevar="time",sep="_",varying=list(c("Obese14","Obese21"),c("Overweight14","Overweight21"),c("hibp14","hibp21")),v.names=c("Obese","Overweight","hibp"),direction="long")
>
> library(car)
> BP.stack3$time<-recode(BP.stack3$time,"1=14;2=21")
> head(BP.stack3,2)
>   #  CODEA Sex MaternalAge Education Birthplace AggScore IntScore  Categ
> time
> #8.1 8   2   4 4  100 Normal
> 14
> #9.1 9   1   3 6  200 Normal
> 14
>   #  Obese Overweight hibp
> #8.1 0  00
>
> Now, your formula: (HIBP~time*Age), is it MaternalAge?
> If it is, it has three values
> unique(BP.stack3$MaternalAge)
> #[1] 4 3 5
> and for time (14,21) # If it says that geese/geeglm, contrasts could be
> applied with factors>=2 levels, what is the problem?
> If you take "Categ" variable, it also has 3 levels (Normal, Obese,
> Overweight).
>
>  BP.stack3$MaternalAge<-factor(BP.stack3$MaternalAge)
>  BP.stack3$time<-factor(BP.stack3$time)
>
> library(geepack)
> For your last question about how to get the p-values:
> # Using one of the example datasets:
> data(seizure)
>  seiz.l <- reshape(seizure,
>varying=list(c("base","y1", "y2", "y3", "y4")),
>v.names="y", times=0:4, direction="long")
>  seiz.l <- seiz.l[order(seiz.l$id, seiz.l$time),]
>  seiz.l$t <- ifelse(seiz.l$time == 0, 8, 2)
>  seiz.l$x <- ifelse(seiz.l$time == 0, 0, 1)
>  m1 <- geese(y ~ offset(log(t)) + x + trt + x:trt, id = id,
>  data=seiz.l, corstr="exch", family=poisson)
>  summary(m1)
>
>  summary(m1)$mean["p"]
> #p
> #(Intercept) 0.000
> #x   0.3347040
> #trt 0.9011982
> #x:trt   0.6236769
>
>
> #If you need the p-values of the scale
>summary(m1)$scale["p"]
>  #   p
> #(Intercept) 0.0254634
>
> Hope it helps.
>
> A.K.
>
>
>
>
>
>
> - Original Message -
> From: rex2013 <[hidden 
> email]>
>
> To: [hidden email] 
> Cc:
> Sent: Sunday, January 6, 2013 4:55 AM
> Subject: Re: [R] random effects model
>
> Hi A.K
>
> Regarding my question on comparing normal/ obese/overweight with 

Re: [R] round 2.3 to 2.5 and 2.2 to 2.0

2013-01-07 Thread Mat
thx, works perfektly :-)



--
View this message in context: 
http://r.789695.n4.nabble.com/round-2-3-to-2-5-and-2-2-to-2-0-tp4654813p4654829.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] Log scale on y axis of parallel coordinate plot (lattice)

2013-01-07 Thread Roland Seubert

Hello all,

I would like to make a parallel coordinate plot with lattice. The plot 
should have vertical log scale axes, and should in principle look like 
this one (I put less chemical elements in my example below):


http://www.geokem.com/images/scans/epr-and-N_Chile_Ridge.gif

The data I am trying to plot are chemical analyses of rock samples (data 
frame "df"). The data needs to be normalised against a reference sample 
(vector "norm"), to get the actual data to be plotted (data frame 
"df_n"). Here is a simplified example:


> df <- data.frame(La = c(3.0, 2.9, 2.7), Eu = c(0.86, 0.76, 0.66), Lu 
= c(0.07, 0.04, 0.04), row.names = c("sample1", "sample2", "sample3"))

> norm <- c(0.237, 0.0563, 0.0246)
> df_n <- df / norm
> df_n
   LaEuLu
sample1  12.65823  3.628692 0.2953586
sample2  51.50977 13.499112 0.7104796
sample3 109.75610 26.829268 1.6260163

The plot needs the same scale for all axes, so my simple panel function 
would be:


> panel.myplot <- function(..., common.scale) {panel.parallel(..., 
common.scale = TRUE)}


I tried to plot the data with the following command to get vertical axes 
with a log scale:


> parallelplot(~ df_n, panel = panel.myplot, horizontal.axis = FALSE, 
scales = list(y = list(log = 10)))


The problem is that lattice simply ignores the log scale and gives me 
the following warning:


Warning message:
In parallelplot.formula(~df_n, panel = panel.spiderplot, horizontal.axis 
= FALSE,  :

  cannot have log y-scale

I browsed through the source code of the respective lattice source file 
and it seems to me that lattice does not support a log scale for y axes 
in parallel coordinate plots? I would greatly appreciate if anyone had 
an idea how to make this plot with lattice or point me to some kind of 
workaround.


Thanks a lot in advance,

Roland

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Predicting New Data -

2013-01-07 Thread rydood
Uwe Ligges-3 wrote
> Hard to debug given we neither have the formula nor the data.

All I can give is that the data is in the form y, x1, x2, x3,... and the
formula is y ~ x1+x2+x3+...
Maybe it's because I'm not explicitly specifying the base learners.



--
View this message in context: 
http://r.789695.n4.nabble.com/Predicting-New-Data-Error-mboost-tp4654677p4654833.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] Refer to previous row

2013-01-07 Thread Paolo Donatelli
Hi all,

I have a very basic doubt -- but still, I am a newby!

My question is about referring to the previous row: in a sample as the
following...

ID  X1  X2
1   A   12
2   A   6
3   A   10
1   B   17
2   B   19
1   C   22
1   D   13
2   D   19
3   D   21

... I would like to create a dummy variable equal to 1 whenever the
value of ID of the current row is lower or equal than the value of ID
of the previous row -- check the new vector X3 I'd like to obtain:

ID  X1  X2  X3
1   A   12  0
2   A   60
3   A   10  0
1   B   17  1
2   B   19  0
1   C   22  1
1   D   13  1
2   D   19  0
3   D   21  0

I have searched a lot without finding a decent and working solution. I
suppose it is just some basic matter of indexing language, something
like

X3<- as.numeric ( ID[n] <= ID[n-1])

but it is not so simple!


thanks!
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] plot x-axis DateTime NOT evenly spaced

2013-01-07 Thread arun
Hi,
Try this:
dat1<-read.table(text="1 2012-07-01 00:57:54 +0900    156
2 2012-07-01 01:07:41 +0900    587
3 2012-07-01 01:09:31 +0900    110
4 2012-07-01 01:18:42 +0900    551
5 2012-07-01 01:39:01 +0900    1219
6 2012-07-01 01:40:40 +0900  99",sep="",header=FALSE,stringsAsFactors=FALSE)

dat2<-data.frame(date=paste(dat1[,2],dat1[,3],paste0("+",dat1[,4]),sep=" 
"),Interval=dat1[,5])
dat2$date<-as.POSIXct(dat2$date,"%Y-%m-%d %H:%M:%S")
library(xts)
 dat3<-xts(dat2[,-1],order.by=dat2[,1])
plot(dat3)
A.K.




- Original Message -
From: ishi soichi 
To: r-help 
Cc: 
Sent: Monday, January 7, 2013 3:55 AM
Subject: [R] plot x-axis DateTime NOT evenly spaced

R-64 latest

Hi. I am trying to plot a set of csv data, which looks like

> head(interval)
                       date inteval
1 2012-07-01 00:57:54 +0900     156
2 2012-07-01 01:07:41 +0900     587
3 2012-07-01 01:09:31 +0900     110
4 2012-07-01 01:18:42 +0900     551
5 2012-07-01 01:39:01 +0900    1219
6 2012-07-01 01:40:40 +0900      99

as you can see, more than one event happens each day, and they are not
evenly spaced.  Obviously hours, minutes and seconds are important for the
plot.

I tried

interval$date <- as.Date(interval$date, "%Y-%m-%d %H:%M:%S +0900")

but this chops the time off.

Could anyone show me how to plot data with x values as Date(or Time)
objects?

soichi

    [[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] posting a question in the R-help forum

2013-01-07 Thread Violet Swakman
Hello,
I wanted to post this question below, on the R-help forum, but I'm not sure
I succeeded because it said that I wasn't subscribed to the mailing list
yet.
Now I am subscribed, but will my question be accepted now automatically, or
should I submit it again?
Thanks in advance,
Violet Swakman



Hello everyone,

I'm having trouble understanding my output from a linear mixed effects
model (nlme :: lme), I hope someone can help me.
Say I'm interested in the effect of Tarsus length on the Bar length of
feathers.
I used an lme since some birds were living in the same territory, so
territory was included as random effect.
Both Bar length and Tarsus length are seen as numerical values, Territory
is seen as a factor.

#
m1 <- lme(Bar_length~Tarsus_length, random = ~ 1|Territory, data=data)

> summary(m1)
Linear mixed-effects model fit by REML
 Data: min_s12
AIC   BIC   logLik
  -104.1593 -99.98116 56.07963

Random effects:
 Formula: ~1 | as.factor(Territory)
(Intercept)   Residual
StdDev:  0.01023884 0.01072872

Fixed effects: Av_bar_length ~ Tarsus_av
  Value  Std.Error   DFt-value p-value
(Intercept)  0.22391092 0.08472658 19  2.6427470  0.0160
Tarsus_av   -0.00048219 0.00338510  2 -0.1424453  0.8998
 Correlation:
  (Intr)
Tarsus_av -0.999

Standardized Within-Group Residuals:
Min Q1Med
Q3 Max
-2.02920116 -0.49093095  0.05736504  0.47632005  1.15944871

Number of Observations: 23
Number of Groups: 20

##

I do not understand why the model needs 17 degrees of freedom to calculate
1 intercept and slope (just 1 numerical explanatory variable). Could anyone
maybe explain this to me?

When I use Season, a factor with 2 levels, as an explanatory variable the
same thing happens, the model takes 17 DF's to calculate the effect of
Season.

Thanks in advance,
Violet

[[alternative HTML version deleted]]

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


Re: [R] how to aggregate T-test result in an elegant way?

2013-01-07 Thread arun
Hi Yao,

It's okay.

How did you generate the 3 D array?
Using ?acast()

I am not sure I understand your question "

if you meet a t-test task as I described  , is that generate a
high-dimension array  a good way ?"

Do you want to do the t-test in the melt dataset?

b<- read.table(text="
ID    O2    variable    value
1    TWF2H5    13% EW.INCU    49.38
2    TWF2H6    13% EW.INCU    48.02
3    TWF2H19    13%    EW.INCU    51.44
280    TWF2H101    13% EW.17.5    42.26
281    TWF2H105    13%    EW.17.5    43.52
282    TWF2H106    13% EW.17.5    42.83
472    TWF2N102    21% EW.17.5    45.97
473    TWF2N104    21%    EW.17.5 43.32
474    TWF2N106    21% EW.17.5    48.63
689    TWF2N2    21% EMW    19.57
690    TWF2N6    21%    EMW    18.07
691    TWF2N10    21%    EMW    15.4
491    TWF2H5    13%    EMW    15.61
492    TWF2H6    13%    EMW    13.41
493    TWF2H19    13%    EMW    14.03
199    TWF2N2    21%    EW.INCU    48.69
200    TWF2N6    21%    EW.INCU    50.52
201    TWF2N10    21%    EW.INCU    42.04
",sep="",header=TRUE,stringsAsFactors=FALSE)
 res<-lapply(lapply(split(b,b$variable),function(x) 
t.test(x$value[x$O2=="13%"],x$value[x$O2=="21%"])),function(x) 
data.frame(mean=x$estimate,p.value=x$p.value))
res1<-do.call(rbind,res)
row.names(res1)[grep("mean of 
x",row.names(res1))]<-gsub("(.*\\.).*$","\\113%",row.names(res1)[grep("mean of 
x",row.names(res1))])
row.names(res1)[grep("mean of 
y",row.names(res1))]<-gsub("(.*\\.).*$","\\121%",row.names(res1)[grep("mean of 
y",row.names(res1))])
res1
#    mean    p.value
#EMW.13% 14.35000 0.09355374
#EMW.21% 17.68000 0.09355374
#EW.17.5.13% 42.87000 0.17464018
#EW.17.5.21% 45.97333 0.17464018
#EW.INCU.13% 49.61333 0.43689727
#EW.INCU.21% 47.08333 0.43689727

A.K.



- Original Message -
From: Yao He 
To: arun 
Cc: R help 
Sent: Monday, January 7, 2013 4:00 AM
Subject: Re: [R] how to aggregate T-test result in an elegant way?

Hi, arun
I'm so sorry for that isn't helpful.
One of question is that I don't know how  to subset a small part as it
is a 3-dimension array so I just show the structure of that.
I tried  dput()  to a file , then what should I do for subsetting it?

Another question is :
My rawdata is a "melt" dataframe like that:
IID    O2    variable    value
1    TWF2H5    13%     EW.INCU    49.38
2    TWF2H6    13%     EW.INCU    48.02
3    TWF2H19    13%     EW.INCU    51.44
280    TWF2H101    13%     EW.17.5    42.26
281    TWF2H105    13%     EW.17.5     43.52
282    TWF2H106    13%     EW.17.5    42.83
472    TWF2N102    21%     EW.17.5    45.97
473    TWF2N104    21%     EW.17.5     43.32
474    TWF2N106    21%     EW.17.5    48.63
689    TWF2N2    21%      EMW    19.57
690    TWF2N6    21%     EMW    18.07
691    TWF2N10    21%    EMW    15.4
491    TWF2H5         13%    EMW    15.61
492    TWF2H6        13%    EMW    13.41
493    TWF2H19    13%    EMW    14.03
199    TWF2N2    21%    EW.INCU    48.69
200    TWF2N6    21%    EW.INCU    50.52
201    TWF2N10    21%    EW.INCU    42.04

if you meet a t-test task as I described  , is that generate a
high-dimension array  a good way ?
Thank you!

Yao He
2013/1/7 arun :
> HI,
> I tried to create an example dataset (as you didn't provide the data).
> set.seed(25)
> a<-array(sample(1:50,60,replace=TRUE),dim=c(2,10,3))
> dimnames(a)[[1]]<-c("13%","21%")
> dimnames(a)[[2]]<-paste("TWF2H",101:110,sep="")
> dimnames(a)[[3]]<-c("EW.INCU","EW.17.5","EMW")
>
>
> str(a)
> # int [1:2, 1:10, 1:3] 21 35 8 45 7 50 32 17 4 15 ...
>  #- attr(*, "dimnames")=List of 3
>   #..$ : chr [1:2] "13%" "21%"
>   .#.$ : chr [1:10] "TWF2H101" "TWF2H102" "TWF2H103" "TWF2H104" ...
>   #..$ : chr [1:3] "EW.INCU" "EW.17.5" "EMW"
>
> res<-lapply(lapply(seq_len(dim(a)[3]),function(i) 
> t.test(a[dimnames(a)[[1]][1],,i],a[dimnames(a)[[1]][2],,i])),function(x) 
> data.frame(mean=x$estimate,p.value=x$p.value))
> res1<-do.call(rbind,res)
>   row.names(res1)[grep("mean of 
>x",row.names(res1))]<-gsub("(.*\\.).*$","\\113%",row.names(res1)[grep("mean of 
>x",row.names(res1))])
>  row.names(res1)[grep("mean of 
>y",row.names(res1))]<-gsub("(.*\\.).*$","\\121%",row.names(res1)[grep("mean of 
>y",row.names(res1))])
> res1
> #            mean   p.value
> #EW.INCU.13% 22.3 0.2754842
> #EW.INCU.21% 29.3 0.2754842
> #EW.17.5.13% 20.5 0.4705772
> #EW.17.5.21% 16.0 0.4705772
> #EMW.13%     23.9 0.9638679
> #EMW.21%     24.2 0.9638679
> A.K.
>
>
>
>
> - Original Message -
> From: Yao He 
> To: arun 
> Cc: R help 
> Sent: Sunday, January 6, 2013 11:21 PM
> Subject: Re: [R] how to aggregate T-test result in an elegant way?
>
> Thank you,it is really helpful everytime.
>
> I didn't provide any example data because I thought it is just a
> question of how to report t.test() result in R.
> However,as you say,it is better to show more details for finding an elegant 
> way
>
> In fact  I generate a 3-dimension array like that:
> str(a)
> num [1:2, 1:245, 1:3] 47.5 NA 48.9 NA 47.5 ..

Re: [R] Predicting New Data -

2013-01-07 Thread Benjamin Hofner
Sorry, but your current information doesn't help to solve your problem. Not
having explicit base-learners is not the cause of your problem, especially
as for glmboost() there are no explicit base-learners. You should definitely
provide a minimal example that helps to reproduce your error/problem. 

For me, everything works as intended:

> ## load mboost
> library(mboost)
> 
> ## set up model (in analogy to ?glmboost)
> cars.gb <- glmboost(dist ~ speed, data = cars,
+ control = boost_control(mstop = 400))
> ## without new data
> p1 <- predict(cars.gb)
> head(p1)
  [,1]
[1,] -1.849460
[2,] -1.849460
[3,]  9.947766
[4,]  9.947766
[5,] 13.880175
[6,] 17.812584
> 
> ## now with new data
> p2 <- predict(cars.gb, newdata =  cars)
> head(p2)
  [,1]
[1,] -1.849460
[2,] -1.849460
[3,]  9.947766
[4,]  9.947766
[5,] 13.880175
[6,] 17.812584
> 
> ## now with really new data
> data <- data.frame(speed = c(1:10))
> p3 <- predict(cars.gb, newdata =  data)
> p3
[,1]
 [1,] -13.646686
 [2,]  -9.714277
 [3,]  -5.781869
 [4,]  -1.849460
 [5,]   2.082949
 [6,]   6.015358
 [7,]   9.947766
 [8,]  13.880175
 [9,]  17.812584
[10,]  21.744993


rydood wrote
> 
> Uwe Ligges-3 wrote
>> Hard to debug given we neither have the formula nor the data.
> All I can give is that the data is in the form y, x1, x2, x3,... and the
> formula is y ~ x1+x2+x3+...
> Maybe it's because I'm not explicitly specifying the base learners.





--
View this message in context: 
http://r.789695.n4.nabble.com/Predicting-New-Data-Error-mboost-tp4654677p4654840.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] R-help post Bayesian CART

2013-01-07 Thread Rachel Plewes
Hi,

I have explored many of the R packages that construct Bayesian trees including 
the tgp, bart, BMA and maptree packages. I have also searched through some 
other packages but they do not seem to be suitable for the type of analysis I 
need to do. I need to construct Bayesian CART that have terminal nodes which 
have bivariate regressions (not multiple regressions like most of the packages 
do). Does anyone know of a package that has this capability?

Thanks,
Rachel
M.Sc Student

[[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] try()-function does not catch error in BATCH-job if Matrix is loaded

2013-01-07 Thread Sarah Brockhaus

Hello,

In my simulation I use the try()-function to catch possible errors when 
fitting models. I run the simulationon a Linux-server using the command 
" R CMD BATCH nameOfFile.R &".  When executing the code as batch-job I 
get the problem that the execution is halted without giving an error 
message. But when I run the code interactivly the error is catched by 
try() as it would be expected.


The problem is somewhat strange as it only occurs when the code is 
executed as a batch-job and when the package "Matrix" is loaded.


I wrote a small example reproducing the error. (In my code the error 
occurs in mgcv:::fixDependence, which looks like the code I'm using 
below in order to get a small reproducible example. I realized that the 
code  makes no sense...)


##
library(Matrix)

R <- matrix(abs(rnorm(25)), 5, 5)
r0 <- r <- nrow(R)

# while-loop produces error that should be catched by the function try()
try(
while (mean(R[r0:r, r0:r]) > 0) r0 <- r0 - 1
)

# so "Hello" should be printed
print("Hello")
##

I use R version 2.15.1 (2012-06-22) on a x86_64-pc-linux-gnu (64-bit) 
platform.


I would be grateful for help.

Best regards,
Sarah Brockhaus

PHD-student
Department of Statistics
University of Munich

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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 aggregate T-test result in an elegant way?

2013-01-07 Thread arun
Hi Yao,

You could also have the results in a wide format:
res<-do.call(rbind,lapply(lapply(split(b,b$variable),function(x) 
t.test(x$value[x$O2=="13%"],x$value[x$O2=="21%"])),function(x) 
data.frame(mean13=x$estimate[1],mean21=x$estimate[2],p.value=x$p.value,CILow=x$conf.int[1],CIHigh=x$conf.int[2])))
 res
#  mean13   mean21    p.value CILow    CIHigh
#EMW 14.35000 17.68000 0.09355374 -7.682686  1.022686
#EW.17.5 42.87000 45.97333 0.17464018 -9.265622  3.058955
#EW.INCU 49.61333 47.08333 0.43689727 -7.119234 12.179234
A.K.




- Original Message -
From: Yao He 
To: arun 
Cc: R help 
Sent: Monday, January 7, 2013 10:57 AM
Subject: Re: [R] how to aggregate T-test result in an elegant way?

Hi,arun

Yes , I just want to do the t.test
I think maybe  it is not necessary to generate a 3D array from the raw
data.frame by acast() at first

Thanks a lot

2013/1/7 arun :
> Hi Yao,
>
> It's okay.
>
> How did you generate the 3 D array?
> Using ?acast()
>
> I am not sure I understand your question "
>
> if you meet a t-test task as I described  , is that generate a
> high-dimension array  a good way ?"
>
> Do you want to do the t-test in the melt dataset?
>
> b<- read.table(text="
> ID    O2    variable    value
> 1    TWF2H5    13%     EW.INCU    49.38
> 2    TWF2H6    13%     EW.INCU    48.02
> 3    TWF2H19    13%    EW.INCU    51.44
> 280    TWF2H101    13%     EW.17.5    42.26
> 281    TWF2H105    13%    EW.17.5    43.52
> 282    TWF2H106    13%     EW.17.5    42.83
> 472    TWF2N102    21%     EW.17.5    45.97
> 473    TWF2N104    21%    EW.17.5     43.32
> 474    TWF2N106    21%     EW.17.5    48.63
> 689    TWF2N2    21%     EMW    19.57
> 690    TWF2N6    21%    EMW    18.07
> 691    TWF2N10    21%    EMW    15.4
> 491    TWF2H5        13%    EMW    15.61
> 492    TWF2H6        13%    EMW    13.41
> 493    TWF2H19    13%    EMW    14.03
> 199    TWF2N2    21%    EW.INCU    48.69
> 200    TWF2N6    21%    EW.INCU    50.52
> 201    TWF2N10    21%    EW.INCU    42.04
> ",sep="",header=TRUE,stringsAsFactors=FALSE)
>  res<-lapply(lapply(split(b,b$variable),function(x) 
>t.test(x$value[x$O2=="13%"],x$value[x$O2=="21%"])),function(x) 
>data.frame(mean=x$estimate,p.value=x$p.value))
> res1<-do.call(rbind,res)
> row.names(res1)[grep("mean of 
> x",row.names(res1))]<-gsub("(.*\\.).*$","\\113%",row.names(res1)[grep("mean 
> of x",row.names(res1))])
> row.names(res1)[grep("mean of 
> y",row.names(res1))]<-gsub("(.*\\.).*$","\\121%",row.names(res1)[grep("mean 
> of y",row.names(res1))])
> res1
> #                mean    p.value
> #EMW.13%     14.35000 0.09355374
> #EMW.21%     17.68000 0.09355374
> #EW.17.5.13% 42.87000 0.17464018
> #EW.17.5.21% 45.97333 0.17464018
> #EW.INCU.13% 49.61333 0.43689727
> #EW.INCU.21% 47.08333 0.43689727
>
> A.K.
>
>
>
> - Original Message -
> From: Yao He 
> To: arun 
> Cc: R help 
> Sent: Monday, January 7, 2013 4:00 AM
> Subject: Re: [R] how to aggregate T-test result in an elegant way?
>
> Hi, arun
> I'm so sorry for that isn't helpful.
> One of question is that I don't know how  to subset a small part as it
> is a 3-dimension array so I just show the structure of that.
> I tried  dput()  to a file , then what should I do for subsetting it?
>
> Another question is :
> My rawdata is a "melt" dataframe like that:
> IID    O2    variable    value
> 1    TWF2H5    13%     EW.INCU    49.38
> 2    TWF2H6    13%     EW.INCU    48.02
> 3    TWF2H19    13%     EW.INCU    51.44
> 280    TWF2H101    13%     EW.17.5    42.26
> 281    TWF2H105    13%     EW.17.5     43.52
> 282    TWF2H106    13%     EW.17.5    42.83
> 472    TWF2N102    21%     EW.17.5    45.97
> 473    TWF2N104    21%     EW.17.5     43.32
> 474    TWF2N106    21%     EW.17.5    48.63
> 689    TWF2N2    21%      EMW    19.57
> 690    TWF2N6    21%     EMW    18.07
> 691    TWF2N10    21%    EMW    15.4
> 491    TWF2H5         13%    EMW    15.61
> 492    TWF2H6        13%    EMW    13.41
> 493    TWF2H19    13%    EMW    14.03
> 199    TWF2N2    21%    EW.INCU    48.69
> 200    TWF2N6    21%    EW.INCU    50.52
> 201    TWF2N10    21%    EW.INCU    42.04
>
> if you meet a t-test task as I described  , is that generate a
> high-dimension array  a good way ?
> Thank you!
>
> Yao He
> 2013/1/7 arun :
>> HI,
>> I tried to create an example dataset (as you didn't provide the data).
>> set.seed(25)
>> a<-array(sample(1:50,60,replace=TRUE),dim=c(2,10,3))
>> dimnames(a)[[1]]<-c("13%","21%")
>> dimnames(a)[[2]]<-paste("TWF2H",101:110,sep="")
>> dimnames(a)[[3]]<-c("EW.INCU","EW.17.5","EMW")
>>
>>
>> str(a)
>> # int [1:2, 1:10, 1:3] 21 35 8 45 7 50 32 17 4 15 ...
>>  #- attr(*, "dimnames")=List of 3
>>   #..$ : chr [1:2] "13%" "21%"
>>   .#.$ : chr [1:10] "TWF2H101" "TWF2H102" "TWF2H103" "TWF2H104" ...
>>   #..$ : chr [1:3] "EW.INCU" "EW.17.5" "EMW"
>>
>> res<-lapply(lapply(seq_len(dim(a)[3]),function(i) 
>> t.test(a[dimnames(a)[[1]][1],,i],a[dimnames(a)[[1]][2],,i])),function(x) 
>> data.frame(mean=x$estimate,

Re: [R] Changing mtext direction, or using text for the margin?

2013-01-07 Thread Michael Rennie


Any thoughts on what that dirty hack might be or any leads on where to 
start? Perhaps a whole new plot region in the margin or something? Is 
that even possible? I'm having a difficult time imagining how I can do this.


Mike

Uwe Ligges wrote:



On 07.01.2013 07:00, Michael Rennie wrote:


Hi all,

I have read through the archives, but can't find a solution to this
problem.

I need the text direction on "dependent B", plotted in margin 4, to go
top to bottom (opposite what it is now). Here's some sample code:

#plot with mtext example

par(mgp = c(2,1,0), mfrow=c(2,2), las=1, mar=c(2,2,2,2), omi=
c(0.5,0.2,0,0.2))

a<-1:10
b<-7:16
c<-21:30

plot(a~b, ylab="", xlab="", xaxt="n")
axis(1, at=b, labels=FALSE, tick=TRUE)
plot(c~b, ylab="", xlab="", xaxt="n", yaxt="n")
axis(1, at=b, labels=FALSE, tick=TRUE)
axis(4, at=c(22,24,26,28,30), labels=TRUE, tick=TRUE)
plot(a~b, ylab="")
plot(c~b, ylab="", yaxt="n")
axis(4, at=c(22,24,26,28,30), labels=TRUE, tick=TRUE)

mtext(side=1, "independent", outer=TRUE, line=1, padj=1)
mtext(side=2, "dependent A", las=0, outer=TRUE, line=0.25)
mtext(side=4, "dependent B", las=0, outer=TRUE, line=0.25)

I have seen an example in help where i can use text and (srt) to
manipulate this with a single panel plot, but not for a multi-panel
example. Any help would be greatly appreciated.



You cannot do that with mtext, you rather need a dirty hack with 
text(), I believe.


Best,
Uwe Ligges




Cheers,

Mike

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] renumber a list of numbers

2013-01-07 Thread Charles Determan Jr
Thank you for your response, I didn't do the -7 because this is just a
small part of the dataset and the number are not all consistent (ie skip
some numbers occasionally).  However, I did not think of the
as.numeric(as.factor()) sequence.  That did the trick, simple lapse in my
thoughts that cost me a lot of time.

Thanks again,
Charles

On Mon, Jan 7, 2013 at 11:23 AM, Sarah Goslee wrote:

> Hi,
>
> It isn't entirely clear what you want, because it seems too simple.
> And most of your sample data are irrelevant, aren't they?
>
> Why not just use:
>
> testdata$Event_name2 <- testdata$Event_name - 7
>
> Or you could try:
>
> testdata$Event_name3 <- as.numeric(as.factor(testdata$Event_name))
>
> which will make the values into consecutive integers.
>Event_name Event_name2 Event_name3
> 9   8   1   1
> 10  9   2   2
> 11 10   3   3
> 12 11   4   4
> 13 12   5   5
> 14 13   6   6
> 15 14   7   7
> 16 15   8   8
> 17 16   9   9
> 18 17  10  10
> 19 18  11  11
> 20 19  12  12
> 21 20  13  13
> 22 21  14  14
> 23 22  15  15
> 24 23  16  16
> 25 24  17  17
> 26 25  18  18
> 27 26  19  19
> 28 27  20  20
> 29 28  21  21
> 30 29  22  22
> 31 30  23  23
> 32 31  24  24
> 33 33  26  25
> 34 34  27  26
> 41  8   1   1
> 42  9   2   2
> 43 10   3   3
> 44 11   4   4
> 45 12   5   5
> 46 13   6   6
> 47 14   7   7
> 48 15   8   8
> 49 16   9   9
> 50 17  10  10
> 51 18  11  11
> 52 19  12  12
> 53 20  13  13
> 54 21  14  14
> 55 22  15  15
> 56 23  16  16
> 57 24  17  17
> 58 25  18  18
> 59 26  19  19
> 60 27  20  20
> 61 28  21  21
> 62 29  22  22
> 63 30  23  23
>
>
> On Mon, Jan 7, 2013 at 11:41 AM, Charles Determan Jr 
> wrote:
> > Greetings R users,
> >
> > I am trying to renumber my groups within the file shown below.  The
> groups
> > are currently set as 8,9,10,etc.  I would like to renumber this as
> > 1,2,3,etc.  I have searched the help files and only come across using the
> > rownames to renumber the values but I need to match values.  Any
> assistance
> > is always appreciated,
> >
> > Regards,
> > Charles
> >
> > structure(list(Group = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L,
> > 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L,
> > 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
> > 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L), .Label = c("NO", "YES"
> > ), class = "factor"), Event_name = c(8, 9, 10, 11, 12, 13, 14,
> > 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30,
> > 31, 33, 34, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20,
> > 21, 22, 23, 24, 25, 26, 27, 28, 29, 30), Glucose.n = c(26, 3,
> > 0, 26, 1, 25, 26, 25, 25, 26, 26, 25, 23, 23, 24, 24, 26, 26,
> > 25, 25, 26, 26, 24, 21, 7, 12, 4, 0, 0, 4, 0, 4, 4, 4, 4, 4,
> > 4, 4, 3, 4, 4, 4, 3, 3, 2, 2, 2, 1, 1), Glucose.m = c(92.5,
> > 90.3,
> > NaN, 97.2307692307692, 116, 97.84, 107.653846153846, 105.32,
> > 102.6, 94.6538461538462, 96.076923076923, 92.24, 87.5652173913043,
> > 79.3913043478261, 81.29167, 77.5, 75.9230769230769,
> > 74.4615384615385,
> > 72.68, 76.32, 74.9615384615385, 72.2307692307692, 92.54167,
> > 105.619047619048, 93.4285714285714, 96.5, 90, NaN, NaN, 86.5,
> > NaN, 87, 90.25, 92.5, 98.75, 95.75, 94, 88.25, 54.3,
> > 52, 74.5, 77.75, 81, 97.3, 82.5, 85, 66.5, 51, 81
> > ), Glucose.sd = c(18.9256439784753, 27.5922694487665, NA,
> 25.3050314242961,
> > NA, 17.3917605012642, 21.027491163127, 12.0094407308029,
> 28.0728219695373,
> > 17.7334538264655, 10.7700439253443, 12.7778454104490, 11.0075432274935,
> > 14.6992242542214, 12.2739709270814, 10.9266328000819, 10.4457573279225,
> > 13.1338669682033, 8.2194890352138, 19.9556174213344, 17.6079090620795,
> > 10.9299869800753, 19.3052801217623, 29.7883806046522, 17.2032665779607,
> > 18.4563366797521, 18.0554700852678, NA, NA, 19.7399763

Re: [R] renumber a list of numbers

2013-01-07 Thread Sarah Goslee
Hi,

It isn't entirely clear what you want, because it seems too simple.
And most of your sample data are irrelevant, aren't they?

Why not just use:

testdata$Event_name2 <- testdata$Event_name - 7

Or you could try:

testdata$Event_name3 <- as.numeric(as.factor(testdata$Event_name))

which will make the values into consecutive integers.
   Event_name Event_name2 Event_name3
9   8   1   1
10  9   2   2
11 10   3   3
12 11   4   4
13 12   5   5
14 13   6   6
15 14   7   7
16 15   8   8
17 16   9   9
18 17  10  10
19 18  11  11
20 19  12  12
21 20  13  13
22 21  14  14
23 22  15  15
24 23  16  16
25 24  17  17
26 25  18  18
27 26  19  19
28 27  20  20
29 28  21  21
30 29  22  22
31 30  23  23
32 31  24  24
33 33  26  25
34 34  27  26
41  8   1   1
42  9   2   2
43 10   3   3
44 11   4   4
45 12   5   5
46 13   6   6
47 14   7   7
48 15   8   8
49 16   9   9
50 17  10  10
51 18  11  11
52 19  12  12
53 20  13  13
54 21  14  14
55 22  15  15
56 23  16  16
57 24  17  17
58 25  18  18
59 26  19  19
60 27  20  20
61 28  21  21
62 29  22  22
63 30  23  23


On Mon, Jan 7, 2013 at 11:41 AM, Charles Determan Jr  wrote:
> Greetings R users,
>
> I am trying to renumber my groups within the file shown below.  The groups
> are currently set as 8,9,10,etc.  I would like to renumber this as
> 1,2,3,etc.  I have searched the help files and only come across using the
> rownames to renumber the values but I need to match values.  Any assistance
> is always appreciated,
>
> Regards,
> Charles
>
> structure(list(Group = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L,
> 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L,
> 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
> 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L), .Label = c("NO", "YES"
> ), class = "factor"), Event_name = c(8, 9, 10, 11, 12, 13, 14,
> 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30,
> 31, 33, 34, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20,
> 21, 22, 23, 24, 25, 26, 27, 28, 29, 30), Glucose.n = c(26, 3,
> 0, 26, 1, 25, 26, 25, 25, 26, 26, 25, 23, 23, 24, 24, 26, 26,
> 25, 25, 26, 26, 24, 21, 7, 12, 4, 0, 0, 4, 0, 4, 4, 4, 4, 4,
> 4, 4, 3, 4, 4, 4, 3, 3, 2, 2, 2, 1, 1), Glucose.m = c(92.5,
> 90.3,
> NaN, 97.2307692307692, 116, 97.84, 107.653846153846, 105.32,
> 102.6, 94.6538461538462, 96.076923076923, 92.24, 87.5652173913043,
> 79.3913043478261, 81.29167, 77.5, 75.9230769230769,
> 74.4615384615385,
> 72.68, 76.32, 74.9615384615385, 72.2307692307692, 92.54167,
> 105.619047619048, 93.4285714285714, 96.5, 90, NaN, NaN, 86.5,
> NaN, 87, 90.25, 92.5, 98.75, 95.75, 94, 88.25, 54.3,
> 52, 74.5, 77.75, 81, 97.3, 82.5, 85, 66.5, 51, 81
> ), Glucose.sd = c(18.9256439784753, 27.5922694487665, NA, 25.3050314242961,
> NA, 17.3917605012642, 21.027491163127, 12.0094407308029, 28.0728219695373,
> 17.7334538264655, 10.7700439253443, 12.7778454104490, 11.0075432274935,
> 14.6992242542214, 12.2739709270814, 10.9266328000819, 10.4457573279225,
> 13.1338669682033, 8.2194890352138, 19.9556174213344, 17.6079090620795,
> 10.9299869800753, 19.3052801217623, 29.7883806046522, 17.2032665779607,
> 18.4563366797521, 18.0554700852678, NA, NA, 19.7399763593239,
> NA, 22.3159136044214, 28.5116935075184, 21.7638844572072, 12.2848144742469,
> 12.0933866224478, 15.5777619273972, 11.842719282327, 39.1066916694999,
> 32.0936130717624, 66.8755062286136, 53.7796429887741, 17.6918060129541,
> 11.5902257671425, 34.6482322781408, 9.89949493661167, 26.1629509039023,
> NA, NA), Glucose.se = c(3.71162415205983, 15.9304041937980, NA,
> 4.96272496248326, NA, 3.47835210025284, 4.12383029856479, 2.40188814616057,
> 5.61456439390745, 3.47781642709786, 2.11217938990701, 2.6908208981,
> 2.29523142628873, 3.06500013246251, 2.50541382409208, 2.23038958057991,
> 2.048

Re: [R] Refer to previous row

2013-01-07 Thread Rui Barradas

Hello,

Try the following.

dat$X3 <- c(0L, dat$ID[-1]  <= dat$ID[-nrow(dat)])

Hope this helps,

Rui Barradas
Em 07-01-2013 13:33, Paolo Donatelli escreveu:

Hi all,

I have a very basic doubt -- but still, I am a newby!

My question is about referring to the previous row: in a sample as the
following...

ID  X1  X2
1   A   12
2   A   6
3   A   10
1   B   17
2   B   19
1   C   22
1   D   13
2   D   19
3   D   21

... I would like to create a dummy variable equal to 1 whenever the
value of ID of the current row is lower or equal than the value of ID
of the previous row -- check the new vector X3 I'd like to obtain:

ID  X1  X2  X3
1   A   12  0
2   A   60
3   A   10  0
1   B   17  1
2   B   19  0
1   C   22  1
1   D   13  1
2   D   19  0
3   D   21  0

I have searched a lot without finding a decent and working solution. I
suppose it is just some basic matter of indexing language, something
like

X3<- as.numeric ( ID[n] <= ID[n-1])

but it is not so simple!


thanks!
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.


__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] During startup - Warning messages:

2013-01-07 Thread Prof Brian Ripley

On 07/01/2013 14:50, Yongjie wrote:

Dear Sir/Madam,

After installing R, I saw this warning message in red when I open the R...

During startup - Warning messages:
1: Setting LC_CTYPE failed, using "C"
2: Setting LC_COLLATE failed, using "C"
3: Setting LC_TIME failed, using "C"
4: Setting LC_MESSAGES failed, using "C"
5: Setting LC_PAPER failed, using "C"

how shall I solve this ?

Thank you very much in advance.


The screenshot is of R.app.  So do follow the posting guide (see the 
footer):


'Platform-specific questions: There are lists R-sig-Mac, R-sig-Debian 
and R-sig-Fedora for R on Mac OS X, Debian/Ubuntu and Fedora/Redhat 
respectively. Questions specific to those platforms (especially re 
installation and the R.app GUI on Mac OS X) are more likely to get 
informed responses on the appropriate list.'


The answer is to set up your Mac correctly, and how to do so is not an R 
question, and the answers will be unintelligible to non-Mac users (and 
probably many Mac users too).  So ask in the right place.




Yj



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




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

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Refer to previous row

2013-01-07 Thread jim holtman
try this


> x <- read.table(text = 'ID  X1  X2
+ 1   A   12
+ 2   A   6
+ 3   A   10
+ 1   B   17
+ 2   B   19
+ 1   C   22
+ 1   D   13
+ 2   D   19
+ 3   D   21', header = TRUE, as.is = TRUE)
> x$X3 <- c(0, diff(x$ID) <= 0)
> x
  ID X1 X2 X3
1  1  A 12  0
2  2  A  6  0
3  3  A 10  0
4  1  B 17  1
5  2  B 19  0
6  1  C 22  1
7  1  D 13  1
8  2  D 19  0
9  3  D 21  0


On Mon, Jan 7, 2013 at 8:33 AM, Paolo Donatelli
 wrote:
> Hi all,
>
> I have a very basic doubt -- but still, I am a newby!
>
> My question is about referring to the previous row: in a sample as the
> following...
>
> ID  X1  X2
> 1   A   12
> 2   A   6
> 3   A   10
> 1   B   17
> 2   B   19
> 1   C   22
> 1   D   13
> 2   D   19
> 3   D   21
>
> ... I would like to create a dummy variable equal to 1 whenever the
> value of ID of the current row is lower or equal than the value of ID
> of the previous row -- check the new vector X3 I'd like to obtain:
>
> ID  X1  X2  X3
> 1   A   12  0
> 2   A   60
> 3   A   10  0
> 1   B   17  1
> 2   B   19  0
> 1   C   22  1
> 1   D   13  1
> 2   D   19  0
> 3   D   21  0
>
> I have searched a lot without finding a decent and working solution. I
> suppose it is just some basic matter of indexing language, something
> like
>
> X3<- as.numeric ( ID[n] <= ID[n-1])
>
> but it is not so simple!
>
>
> thanks!
> 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.



-- 
Jim Holtman
Data Munger Guru

What is the problem that you are trying to solve?
Tell me what you want to do, not how you want to do it.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] During startup - Warning messages:

2013-01-07 Thread R. Michael Weylandt
I'd suggest you do as the following message suggests and follow the
R-MAC-FAQ section 9 notes on locales and whatnot. If you need specific
advice on following this, we'll need to know what version of OS X
you're running and what your local settings are.

That said, I don't think you'll actually hit too much trouble if you
stick to ASCII characters.

MW

On Mon, Jan 7, 2013 at 2:50 PM, Yongjie  wrote:
> Dear Sir/Madam,
>
> After installing R, I saw this warning message in red when I open the R...
>
> During startup - Warning messages:
> 1: Setting LC_CTYPE failed, using "C"
> 2: Setting LC_COLLATE failed, using "C"
> 3: Setting LC_TIME failed, using "C"
> 4: Setting LC_MESSAGES failed, using "C"
> 5: Setting LC_PAPER failed, using "C"
>
> how shall I solve this ?
>
> Thank you very much in advance.
>
> Yj
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/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] Refer to previous row

2013-01-07 Thread David Winsemius


On Jan 7, 2013, at 5:33 AM, Paolo Donatelli wrote:


Hi all,

I have a very basic doubt -- but still, I am a newby!

My question is about referring to the previous row: in a sample as the
following...

ID  X1  X2
1   A   12
2   A   6
3   A   10
1   B   17
2   B   19
1   C   22
1   D   13
2   D   19
3   D   21

... I would like to create a dummy variable equal to 1 whenever the
value of ID of the current row is lower or equal than the value of ID
of the previous row -- check the new vector X3 I'd like to obtain:

ID  X1  X2  X3
1   A   12  0
2   A   60
3   A   10  0
1   B   17  1
2   B   19  0
1   C   22  1
1   D   13  1
2   D   19  0
3   D   21  0


Something like (untested):

  dfrm$X3 <- c(0, as.numeric( diff(dfrm$ID) <= 0 ) )

That might be faster than this sort of untested strategy:

   ...   <- with(dfrm, c( 0 , as.numeric( ID[2:nrow(dfrm)] <= ID[1: 
(nrow(dfrm)-1] ) ) )


In my newbie days I thought a function named `lag` would do it, but  
discovered it was only working on ts-class objects.





I have searched a lot without finding a decent and working solution.  
I Adding

suppose it is just some basic matter of indexing language, something
like

X3<- as.numeric ( ID[n] <= ID[n-1])


An explicit sequence rather than using mathematical notation is  
needed. And if you are using dataframes, you should not be using  
`attach`. That X3 would not be constructed in the dataframe.


--

David Winsemius, MD
Alameda, 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] Refer to previous row

2013-01-07 Thread Berend Hasselman

On 07-01-2013, at 14:33, Paolo Donatelli  wrote:

> Hi all,
> 
> I have a very basic doubt -- but still, I am a newby!
> 
> My question is about referring to the previous row: in a sample as the
> following...
> 
> ID  X1  X2
> 1   A   12
> 2   A   6
> 3   A   10
> 1   B   17
> 2   B   19
> 1   C   22
> 1   D   13
> 2   D   19
> 3   D   21
> 
> ... I would like to create a dummy variable equal to 1 whenever the
> value of ID of the current row is lower or equal than the value of ID
> of the previous row -- check the new vector X3 I'd like to obtain:
> 
> ID  X1  X2  X3
> 1   A   12  0
> 2   A   60
> 3   A   10  0
> 1   B   17  1
> 2   B   19  0
> 1   C   22  1
> 1   D   13  1
> 2   D   19  0
> 3   D   21  0
> 
> I have searched a lot without finding a decent and working solution. I
> suppose it is just some basic matter of indexing language, something
> like
> 
> X3<- as.numeric ( ID[n] <= ID[n-1])

Assuming your data are in a data.frame called dat this will work

dat$X3 <- c(0, diff(dat$ID)<=0)

Berend

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


Re: [R] Changing mtext direction, or using text for the margin?

2013-01-07 Thread Uwe Ligges



On 07.01.2013 17:06, Michael Rennie wrote:


Any thoughts on what that dirty hack might be or any leads on where to
start? Perhaps a whole new plot region in the margin or something? Is
that even possible? I'm having a difficult time imagining how I can do
this.


Ideas fir ugly hacks:

Either use text() (and you rely on the spacing of the current device) or 
try to use layout() and small figure at the reight hand size of all the 
others.


Or go to the grid package and make it less ugly.

Best,
Uwe Ligges



Mike

Uwe Ligges wrote:



On 07.01.2013 07:00, Michael Rennie wrote:


Hi all,

I have read through the archives, but can't find a solution to this
problem.

I need the text direction on "dependent B", plotted in margin 4, to go
top to bottom (opposite what it is now). Here's some sample code:

#plot with mtext example

par(mgp = c(2,1,0), mfrow=c(2,2), las=1, mar=c(2,2,2,2), omi=
c(0.5,0.2,0,0.2))

a<-1:10
b<-7:16
c<-21:30

plot(a~b, ylab="", xlab="", xaxt="n")
axis(1, at=b, labels=FALSE, tick=TRUE)
plot(c~b, ylab="", xlab="", xaxt="n", yaxt="n")
axis(1, at=b, labels=FALSE, tick=TRUE)
axis(4, at=c(22,24,26,28,30), labels=TRUE, tick=TRUE)
plot(a~b, ylab="")
plot(c~b, ylab="", yaxt="n")
axis(4, at=c(22,24,26,28,30), labels=TRUE, tick=TRUE)

mtext(side=1, "independent", outer=TRUE, line=1, padj=1)
mtext(side=2, "dependent A", las=0, outer=TRUE, line=0.25)
mtext(side=4, "dependent B", las=0, outer=TRUE, line=0.25)

I have seen an example in help where i can use text and (srt) to
manipulate this with a single panel plot, but not for a multi-panel
example. Any help would be greatly appreciated.



You cannot do that with mtext, you rather need a dirty hack with
text(), I believe.

Best,
Uwe Ligges




Cheers,

Mike

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] renumber a list of numbers

2013-01-07 Thread Rui Barradas

Hello,

I'm not sure I understand. Do you want to renumber column Event_name 
starting at 1? If so the following does the job.


dat$Event_name <- dat$Event_name - min(dat$Event_name) + 1


Hope this helps,

Rui Barradas
Em 07-01-2013 16:41, Charles Determan Jr escreveu:

Greetings R users,

I am trying to renumber my groups within the file shown below.  The groups
are currently set as 8,9,10,etc.  I would like to renumber this as
1,2,3,etc.  I have searched the help files and only come across using the
rownames to renumber the values but I need to match values.  Any assistance
is always appreciated,

Regards,
Charles

structure(list(Group = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L,
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L,
1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L), .Label = c("NO", "YES"
), class = "factor"), Event_name = c(8, 9, 10, 11, 12, 13, 14,
15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30,
31, 33, 34, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20,
21, 22, 23, 24, 25, 26, 27, 28, 29, 30), Glucose.n = c(26, 3,
0, 26, 1, 25, 26, 25, 25, 26, 26, 25, 23, 23, 24, 24, 26, 26,
25, 25, 26, 26, 24, 21, 7, 12, 4, 0, 0, 4, 0, 4, 4, 4, 4, 4,
4, 4, 3, 4, 4, 4, 3, 3, 2, 2, 2, 1, 1), Glucose.m = c(92.5,
90.3,
NaN, 97.2307692307692, 116, 97.84, 107.653846153846, 105.32,
102.6, 94.6538461538462, 96.076923076923, 92.24, 87.5652173913043,
79.3913043478261, 81.29167, 77.5, 75.9230769230769,
74.4615384615385,
72.68, 76.32, 74.9615384615385, 72.2307692307692, 92.54167,
105.619047619048, 93.4285714285714, 96.5, 90, NaN, NaN, 86.5,
NaN, 87, 90.25, 92.5, 98.75, 95.75, 94, 88.25, 54.3,
52, 74.5, 77.75, 81, 97.3, 82.5, 85, 66.5, 51, 81
), Glucose.sd = c(18.9256439784753, 27.5922694487665, NA, 25.3050314242961,
NA, 17.3917605012642, 21.027491163127, 12.0094407308029, 28.0728219695373,
17.7334538264655, 10.7700439253443, 12.7778454104490, 11.0075432274935,
14.6992242542214, 12.2739709270814, 10.9266328000819, 10.4457573279225,
13.1338669682033, 8.2194890352138, 19.9556174213344, 17.6079090620795,
10.9299869800753, 19.3052801217623, 29.7883806046522, 17.2032665779607,
18.4563366797521, 18.0554700852678, NA, NA, 19.7399763593239,
NA, 22.3159136044214, 28.5116935075184, 21.7638844572072, 12.2848144742469,
12.0933866224478, 15.5777619273972, 11.842719282327, 39.1066916694999,
32.0936130717624, 66.8755062286136, 53.7796429887741, 17.6918060129541,
11.5902257671425, 34.6482322781408, 9.89949493661167, 26.1629509039023,
NA, NA), Glucose.se = c(3.71162415205983, 15.9304041937980, NA,
4.96272496248326, NA, 3.47835210025284, 4.12383029856479, 2.40188814616057,
5.61456439390745, 3.47781642709786, 2.11217938990701, 2.6908208981,
2.29523142628873, 3.06500013246251, 2.50541382409208, 2.23038958057991,
2.04858155574352, 2.57576322922309, 1.64389780704276, 3.99112348426689,
3.45319507311964, 2.14354680364304, 3.94067380331773, 6.50035756909396,
6.50222358617618, 5.32788547515463, 9.0277350426339, NA, NA,
9.86998817966195, NA, 11.1579568022107, 14.2558467537592, 10.8819422286036,
6.14240723712346, 6.04669331122391, 7.7096369861, 5.92135964116351,
22.5782589625015, 16.0468065358812, 33.4377531143068, 26.8898214943871,
10.2143689640297, 6.69161996662824, 24.5, 7, 18.5, NA, NA)), .Names =
c("Group",
"Event_name", "Glucose.n", "Glucose.m", "Glucose.sd", "Glucose.se"
), row.names = c(9L, 10L, 11L, 12L, 13L, 14L, 15L, 16L, 17L,
18L, 19L, 20L, 21L, 22L, 23L, 24L, 25L, 26L, 27L, 28L, 29L, 30L,
31L, 32L, 33L, 34L, 41L, 42L, 43L, 44L, 45L, 46L, 47L, 48L, 49L,
50L, 51L, 52L, 53L, 54L, 55L, 56L, 57L, 58L, 59L, 60L, 61L, 62L,
63L), class = "data.frame")

[[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] During startup - Warning messages:

2013-01-07 Thread David Winsemius


On Jan 7, 2013, at 6:50 AM, Yongjie wrote:


Dear Sir/Madam,

After installing R, I saw this warning message in red when I open  
the R...


During startup - Warning messages:
1: Setting LC_CTYPE failed, using "C"
2: Setting LC_COLLATE failed, using "C"
3: Setting LC_TIME failed, using "C"
4: Setting LC_MESSAGES failed, using "C"
5: Setting LC_PAPER failed, using "C"

how shall I solve this ?

Thank you very much in advance.

Yj
__


Rather than pasting screenshots you should copy the console output to  
you message body.


Have you read and followed the suggestions in the warning message to  
read the MacOS RFAQ? Its accessible in a GUI session from the "Help"  
menu.


You should at the very least say what your prefered language and  
timezone settings might be.



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


Read the Posting Guide, too.


and provide commented, minimal, self-contained, reproducible code.



--
David Winsemius, MD
Alameda, 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] renumber a list of numbers

2013-01-07 Thread David Winsemius


On Jan 7, 2013, at 8:41 AM, Charles Determan Jr wrote:


Greetings R users,

I am trying to renumber my groups within the file shown below.  The  
groups

are currently set as 8,9,10,etc.  I would like to renumber this as
1,2,3,etc.  I have searched the help files and only come across  
using the
rownames to renumber the values but I need to match values.  Any  
assistance

is always appreciated,


Assuming you mean the `Event_name` column rather than the `Groups`  
column which currently is a factor with label values of "YES" and  
"NO", then why not:


dfrm$Event_name column <- dfrm$Event_name - 7

-- David.



Regards,
Charles

structure(list(Group = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L,
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L,
1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L), .Label = c("NO", "YES"
), class = "factor"), Event_name = c(8, 9, 10, 11, 12, 13, 14,
15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30,
31, 33, 34, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20,
21, 22, 23, 24, 25, 26, 27, 28, 29, 30), Glucose.n = c(26, 3,
0, 26, 1, 25, 26, 25, 25, 26, 26, 25, 23, 23, 24, 24, 26, 26,
25, 25, 26, 26, 24, 21, 7, 12, 4, 0, 0, 4, 0, 4, 4, 4, 4, 4,
4, 4, 3, 4, 4, 4, 3, 3, 2, 2, 2, 1, 1), Glucose.m = c(92.5,
90.3,
NaN, 97.2307692307692, 116, 97.84, 107.653846153846, 105.32,
102.6, 94.6538461538462, 96.076923076923, 92.24, 87.5652173913043,
79.3913043478261, 81.29167, 77.5, 75.9230769230769,
74.4615384615385,
72.68, 76.32, 74.9615384615385, 72.2307692307692, 92.54167,
105.619047619048, 93.4285714285714, 96.5, 90, NaN, NaN, 86.5,
NaN, 87, 90.25, 92.5, 98.75, 95.75, 94, 88.25, 54.3,
52, 74.5, 77.75, 81, 97.3, 82.5, 85, 66.5, 51, 81
), Glucose.sd = c(18.9256439784753, 27.5922694487665, NA,  
25.3050314242961,
NA, 17.3917605012642, 21.027491163127, 12.0094407308029,  
28.0728219695373,
17.7334538264655, 10.7700439253443, 12.7778454104490,  
11.0075432274935,
14.6992242542214, 12.2739709270814, 10.9266328000819,  
10.4457573279225,

13.1338669682033, 8.2194890352138, 19.9556174213344, 17.6079090620795,
10.9299869800753, 19.3052801217623, 29.7883806046522,  
17.2032665779607,

18.4563366797521, 18.0554700852678, NA, NA, 19.7399763593239,
NA, 22.3159136044214, 28.5116935075184, 21.7638844572072,  
12.2848144742469,

12.0933866224478, 15.5777619273972, 11.842719282327, 39.1066916694999,
32.0936130717624, 66.8755062286136, 53.7796429887741,  
17.6918060129541,
11.5902257671425, 34.6482322781408, 9.89949493661167,  
26.1629509039023,

NA, NA), Glucose.se = c(3.71162415205983, 15.9304041937980, NA,
4.96272496248326, NA, 3.47835210025284, 4.12383029856479,  
2.40188814616057,
5.61456439390745, 3.47781642709786, 2.11217938990701,  
2.6908208981,
2.29523142628873, 3.06500013246251, 2.50541382409208,  
2.23038958057991,
2.04858155574352, 2.57576322922309, 1.64389780704276,  
3.99112348426689,
3.45319507311964, 2.14354680364304, 3.94067380331773,  
6.50035756909396,

6.50222358617618, 5.32788547515463, 9.0277350426339, NA, NA,
9.86998817966195, NA, 11.1579568022107, 14.2558467537592,  
10.8819422286036,
6.14240723712346, 6.04669331122391, 7.7096369861,  
5.92135964116351,
22.5782589625015, 16.0468065358812, 33.4377531143068,  
26.8898214943871,

10.2143689640297, 6.69161996662824, 24.5, 7, 18.5, NA, NA)), .Names =
c("Group",
"Event_name", "Glucose.n", "Glucose.m", "Glucose.sd", "Glucose.se"
), row.names = c(9L, 10L, 11L, 12L, 13L, 14L, 15L, 16L, 17L,
18L, 19L, 20L, 21L, 22L, 23L, 24L, 25L, 26L, 27L, 28L, 29L, 30L,
31L, 32L, 33L, 34L, 41L, 42L, 43L, 44L, 45L, 46L, 47L, 48L, 49L,
50L, 51L, 52L, 53L, 54L, 55L, 56L, 57L, 58L, 59L, 60L, 61L, 62L,
63L), class = "data.frame")

[[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
Alameda, 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] Refer to previous row

2013-01-07 Thread Duncan Murdoch

On 07/01/2013 8:33 AM, Paolo Donatelli wrote:

Hi all,

I have a very basic doubt -- but still, I am a newby!

My question is about referring to the previous row: in a sample as the
following...

ID  X1  X2
1   A   12
2   A   6
3   A   10
1   B   17
2   B   19
1   C   22
1   D   13
2   D   19
3   D   21

... I would like to create a dummy variable equal to 1 whenever the
value of ID of the current row is lower or equal than the value of ID
of the previous row -- check the new vector X3 I'd like to obtain:

ID  X1  X2  X3
1   A   12  0
2   A   60
3   A   10  0
1   B   17  1
2   B   19  0
1   C   22  1
1   D   13  1
2   D   19  0
3   D   21  0

I have searched a lot without finding a decent and working solution. I
suppose it is just some basic matter of indexing language, something
like

X3<- as.numeric ( ID[n] <= ID[n-1])

but it is not so simple!


Negative indexing lets you leave out an entry, so x[-1] leaves out the 
first entry, and x[-length(x)] leaves out the last one.  To talk about 
previous entries, you need to do something about the fact that the first 
row has no previous entry.  You gave X3[1] the value 0, suggesting that 
you want to implicitly have the "zeroth" row to have the smallest 
possible value.  So


prevID <- c( -Inf, ID[-length(ID)] )
X3 <- as.numeric( ID <= prevID )

should do what you want.

Duncan Murdoch

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] renumber a list of numbers

2013-01-07 Thread David Winsemius


On Jan 7, 2013, at 10:13 AM, David Winsemius wrote:



On Jan 7, 2013, at 8:41 AM, Charles Determan Jr wrote:


Greetings R users,

I am trying to renumber my groups within the file shown below.  The  
groups

are currently set as 8,9,10,etc.  I would like to renumber this as
1,2,3,etc.  I have searched the help files and only come across  
using the
rownames to renumber the values but I need to match values.  Any  
assistance

is always appreciated,


Assuming you mean the `Event_name` column rather than the `Groups`  
column which currently is a factor with label values of "YES" and  
"NO", then why not:


dfrm$Event_name column <- dfrm$Event_name - 7


Ooops. Didn't mean to leave 'column' in the code. Make that:

dfrm$Event_name <- dfrm$Event_name - 7



-- David.



Regards,
Charles

structure(list(Group = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L,
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L,
1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L), .Label = c("NO", "YES"
), class = "factor"), Event_name = c(8, 9, 10, 11, 12, 13, 14,
15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30,
31, 33, 34, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20,
21, 22, 23, 24, 25, 26, 27, 28, 29, 30), Glucose.n = c(26, 3,
0, 26, 1, 25, 26, 25, 25, 26, 26, 25, 23, 23, 24, 24, 26, 26,
25, 25, 26, 26, 24, 21, 7, 12, 4, 0, 0, 4, 0, 4, 4, 4, 4, 4,
4, 4, 3, 4, 4, 4, 3, 3, 2, 2, 2, 1, 1), Glucose.m = c(92.5,
90.3,
NaN, 97.2307692307692, 116, 97.84, 107.653846153846, 105.32,
102.6, 94.6538461538462, 96.076923076923, 92.24, 87.5652173913043,
79.3913043478261, 81.29167, 77.5, 75.9230769230769,
74.4615384615385,
72.68, 76.32, 74.9615384615385, 72.2307692307692, 92.54167,
105.619047619048, 93.4285714285714, 96.5, 90, NaN, NaN, 86.5,
NaN, 87, 90.25, 92.5, 98.75, 95.75, 94, 88.25, 54.3,
52, 74.5, 77.75, 81, 97.3, 82.5, 85, 66.5, 51, 81
), Glucose.sd = c(18.9256439784753, 27.5922694487665, NA,  
25.3050314242961,
NA, 17.3917605012642, 21.027491163127, 12.0094407308029,  
28.0728219695373,
17.7334538264655, 10.7700439253443, 12.7778454104490,  
11.0075432274935,
14.6992242542214, 12.2739709270814, 10.9266328000819,  
10.4457573279225,
13.1338669682033, 8.2194890352138, 19.9556174213344,  
17.6079090620795,
10.9299869800753, 19.3052801217623, 29.7883806046522,  
17.2032665779607,

18.4563366797521, 18.0554700852678, NA, NA, 19.7399763593239,
NA, 22.3159136044214, 28.5116935075184, 21.7638844572072,  
12.2848144742469,
12.0933866224478, 15.5777619273972, 11.842719282327,  
39.1066916694999,
32.0936130717624, 66.8755062286136, 53.7796429887741,  
17.6918060129541,
11.5902257671425, 34.6482322781408, 9.89949493661167,  
26.1629509039023,

NA, NA), Glucose.se = c(3.71162415205983, 15.9304041937980, NA,
4.96272496248326, NA, 3.47835210025284, 4.12383029856479,  
2.40188814616057,
5.61456439390745, 3.47781642709786, 2.11217938990701,  
2.6908208981,
2.29523142628873, 3.06500013246251, 2.50541382409208,  
2.23038958057991,
2.04858155574352, 2.57576322922309, 1.64389780704276,  
3.99112348426689,
3.45319507311964, 2.14354680364304, 3.94067380331773,  
6.50035756909396,

6.50222358617618, 5.32788547515463, 9.0277350426339, NA, NA,
9.86998817966195, NA, 11.1579568022107, 14.2558467537592,  
10.8819422286036,
6.14240723712346, 6.04669331122391, 7.7096369861,  
5.92135964116351,
22.5782589625015, 16.0468065358812, 33.4377531143068,  
26.8898214943871,

10.2143689640297, 6.69161996662824, 24.5, 7, 18.5, NA, NA)), .Names =
c("Group",
"Event_name", "Glucose.n", "Glucose.m", "Glucose.sd", "Glucose.se"
), row.names = c(9L, 10L, 11L, 12L, 13L, 14L, 15L, 16L, 17L,
18L, 19L, 20L, 21L, 22L, 23L, 24L, 25L, 26L, 27L, 28L, 29L, 30L,
31L, 32L, 33L, 34L, 41L, 42L, 43L, 44L, 45L, 46L, 47L, 48L, 49L,
50L, 51L, 52L, 53L, 54L, 55L, 56L, 57L, 58L, 59L, 60L, 61L, 62L,
63L), class = "data.frame")

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


David Winsemius, MD
Alameda, 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.


[R] R: Re: Refer to previous row

2013-01-07 Thread Paolo Donatelli
It works!

The rationale, if I have understood well, is to take my vector of N elements, 
ask to remove the first/last element, and replace the blank space with a given 
value.

Thank you all for the support!
-Original Message-
From: Duncan Murdoch 
Date: Mon, 07 Jan 2013 13:16:50 
To: 
Cc: 
Subject: Re: [R] Refer to previous row

On 07/01/2013 8:33 AM, Paolo Donatelli wrote:
> Hi all,
>
> I have a very basic doubt -- but still, I am a newby!
>
> My question is about referring to the previous row: in a sample as the
> following...
>
> ID  X1  X2
> 1   A   12
> 2   A   6
> 3   A   10
> 1   B   17
> 2   B   19
> 1   C   22
> 1   D   13
> 2   D   19
> 3   D   21
>
> ... I would like to create a dummy variable equal to 1 whenever the
> value of ID of the current row is lower or equal than the value of ID
> of the previous row -- check the new vector X3 I'd like to obtain:
>
> ID  X1  X2  X3
> 1   A   12  0
> 2   A   60
> 3   A   10  0
> 1   B   17  1
> 2   B   19  0
> 1   C   22  1
> 1   D   13  1
> 2   D   19  0
> 3   D   21  0
>
> I have searched a lot without finding a decent and working solution. I
> suppose it is just some basic matter of indexing language, something
> like
>
> X3<- as.numeric ( ID[n] <= ID[n-1])
>
> but it is not so simple!

Negative indexing lets you leave out an entry, so x[-1] leaves out the 
first entry, and x[-length(x)] leaves out the last one.  To talk about 
previous entries, you need to do something about the fact that the first 
row has no previous entry.  You gave X3[1] the value 0, suggesting that 
you want to implicitly have the "zeroth" row to have the smallest 
possible value.  So

prevID <- c( -Inf, ID[-length(ID)] )
X3 <- as.numeric( ID <= prevID )

should do what you want.

Duncan Murdoch


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


Re: [R] Changing mtext direction, or using text for the margin?

2013-01-07 Thread David Winsemius


On Jan 7, 2013, at 8:06 AM, Michael Rennie wrote:



Any thoughts on what that dirty hack might be or any leads on where  
to start? Perhaps a whole new plot region in the margin or  
something? Is that even possible? I'm having a difficult time  
imagining how I can do this.


The `text` function accepts the 'xpd' argument for `par`.

?par
?text

(Similar but not exactly the same as R-FAQ "7.27 How can I create  
rotated axis labels?")


I'm not sure what sort of "dirt" Uwe was expecting, but I'm guessing  
he meant that you probably were not going to be getting exactly what  
you wanted with the first attempt and were going to need to adjust the  
locations "by hand" after digging around in the par documentation.


--
David.

Mike

Uwe Ligges wrote:



On 07.01.2013 07:00, Michael Rennie wrote:


Hi all,

I have read through the archives, but can't find a solution to this
problem.

I need the text direction on "dependent B", plotted in margin 4,  
to go

top to bottom (opposite what it is now). Here's some sample code:

#plot with mtext example

par(mgp = c(2,1,0), mfrow=c(2,2), las=1, mar=c(2,2,2,2), omi=
c(0.5,0.2,0,0.2))

a<-1:10
b<-7:16
c<-21:30

plot(a~b, ylab="", xlab="", xaxt="n")
axis(1, at=b, labels=FALSE, tick=TRUE)
plot(c~b, ylab="", xlab="", xaxt="n", yaxt="n")
axis(1, at=b, labels=FALSE, tick=TRUE)
axis(4, at=c(22,24,26,28,30), labels=TRUE, tick=TRUE)
plot(a~b, ylab="")
plot(c~b, ylab="", yaxt="n")
axis(4, at=c(22,24,26,28,30), labels=TRUE, tick=TRUE)

mtext(side=1, "independent", outer=TRUE, line=1, padj=1)
mtext(side=2, "dependent A", las=0, outer=TRUE, line=0.25)
mtext(side=4, "dependent B", las=0, outer=TRUE, line=0.25)

I have seen an example in help where i can use text and (srt) to
manipulate this with a single panel plot, but not for a multi-panel
example. Any help would be greatly appreciated.



You cannot do that with mtext, you rather need a dirty hack with  
text(), I believe.


Best,
Uwe Ligges




Cheers,

Mike

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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.


David Winsemius, MD
Alameda, 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] Refer to previous row

2013-01-07 Thread William Dunlap
> I suppose it is just some basic matter of indexing language, something like
> 
> X3<- as.numeric ( ID[n] <= ID[n-1])
> 
> but it is not so simple!

If you first define 'n' as
   n <- seq_along(ID)[-1] # 2:length(ID)
then that code works.

> with(Data, { n <- seq_along(ID)[-1] ; as.numeric ( ID[n] <= ID[n-1]) }) 
[1] 0 0 1 0 1 1 0 0

You may want to put a NA at the start of the result so it has
the same length as the other columns in the data.frame.

I usually find it more convenient to leave things like ID[n]<=ID[n-1]
as logical variables instead of converting them to numbers.  Most
functions do the conversion implicitly if required and leaving them
as logicals helps me remember the meaning.

Bill Dunlap
Spotfire, TIBCO Software
wdunlap tibco.com


> -Original Message-
> From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On 
> Behalf
> Of Paolo Donatelli
> Sent: Monday, January 07, 2013 5:34 AM
> To: r-help@r-project.org
> Subject: [R] Refer to previous row
> 
> Hi all,
> 
> I have a very basic doubt -- but still, I am a newby!
> 
> My question is about referring to the previous row: in a sample as the
> following...
> 
> ID  X1  X2
> 1   A   12
> 2   A   6
> 3   A   10
> 1   B   17
> 2   B   19
> 1   C   22
> 1   D   13
> 2   D   19
> 3   D   21
> 
> ... I would like to create a dummy variable equal to 1 whenever the
> value of ID of the current row is lower or equal than the value of ID
> of the previous row -- check the new vector X3 I'd like to obtain:
> 
> ID  X1  X2  X3
> 1   A   12  0
> 2   A   60
> 3   A   10  0
> 1   B   17  1
> 2   B   19  0
> 1   C   22  1
> 1   D   13  1
> 2   D   19  0
> 3   D   21  0
> 
> I have searched a lot without finding a decent and working solution. I
> suppose it is just some basic matter of indexing language, something
> like
> 
> X3<- as.numeric ( ID[n] <= ID[n-1])
> 
> but it is not so simple!
> 
> 
> thanks!
> 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.

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


Re: [R] Changing mtext direction, or using text for the margin?

2013-01-07 Thread Michael Rennie


Right- I think the issue is getting a text() command (which does accept 
srt) to work outside the plot region, and in the margin where I need it 
to appear. This is easy to "hack" in a single plot, but not so clear how 
to do it in a multi-panel environment with an outer margin.


I suspect going to grid() is where I'm going to have to go (after 
reading the relevant chapters in Paul Murrell's book)


Mike

David Winsemius wrote:


On Jan 7, 2013, at 8:06 AM, Michael Rennie wrote:



Any thoughts on what that dirty hack might be or any leads on where 
to start? Perhaps a whole new plot region in the margin or something? 
Is that even possible? I'm having a difficult time imagining how I 
can do this.


The `text` function accepts the 'xpd' argument for `par`.

?par
?text

(Similar but not exactly the same as R-FAQ "7.27 How can I create 
rotated axis labels?")


I'm not sure what sort of "dirt" Uwe was expecting, but I'm guessing 
he meant that you probably were not going to be getting exactly what 
you wanted with the first attempt and were going to need to adjust the 
locations "by hand" after digging around in the par documentation.




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


Re: [R] Changing mtext direction, or using text for the margin?

2013-01-07 Thread David Winsemius


On Jan 7, 2013, at 10:45 AM, Michael Rennie wrote:



Right- I think the issue is getting a text() command (which does  
accept srt) to work outside the plot region, and in the margin where  
I need it to appear. This is easy to "hack" in a single plot, but  
not so clear how to do it in a multi-panel environment with an outer  
margin.


I suspect going to grid() is where I'm going to have to go (after  
reading the relevant chapters in Paul Murrell's book)




If you are going that route you may want to look at the gridBase  
package.


--
David.

Mike

David Winsemius wrote:


On Jan 7, 2013, at 8:06 AM, Michael Rennie wrote:



Any thoughts on what that dirty hack might be or any leads on  
where to start? Perhaps a whole new plot region in the margin or  
something? Is that even possible? I'm having a difficult time  
imagining how I can do this.


The `text` function accepts the 'xpd' argument for `par`.

?par
?text

(Similar but not exactly the same as R-FAQ "7.27 How can I create  
rotated axis labels?")


I'm not sure what sort of "dirt" Uwe was expecting, but I'm  
guessing he meant that you probably were not going to be getting  
exactly what you wanted with the first attempt and were going to  
need to adjust the locations "by hand" after digging around in the  
par documentation.






David Winsemius, MD
Alameda, 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] Changing mtext direction, or using text for the margin?

2013-01-07 Thread ilai
On Mon, Jan 7, 2013 at 11:52 AM, David Winsemius wrote:

>
> If you are going that route you may want to look at the gridBase package.
>
> Yes for mixing base and grid graphics but IMHO overkill here. Replacing
the last mtext line with

grid::grid.text('dependent B', 0.985 , 0.5 , rot = 270)

should take care of it (leaving exact placement and justification to the OP)

cheers



> --
> David.
>
>  Mike
>>
>> David Winsemius wrote:
>>
>>>
>>> On Jan 7, 2013, at 8:06 AM, Michael Rennie wrote:
>>>
>>>
 Any thoughts on what that dirty hack might be or any leads on where to
 start? Perhaps a whole new plot region in the margin or something? Is that
 even possible? I'm having a difficult time imagining how I can do this.

>>>
>>> The `text` function accepts the 'xpd' argument for `par`.
>>>
>>> ?par
>>> ?text
>>>
>>> (Similar but not exactly the same as R-FAQ "7.27 How can I create
>>> rotated axis labels?")
>>>
>>> I'm not sure what sort of "dirt" Uwe was expecting, but I'm guessing he
>>> meant that you probably were not going to be getting exactly what you
>>> wanted with the first attempt and were going to need to adjust the
>>> locations "by hand" after digging around in the par documentation.
>>>
>>>
>>
> David Winsemius, MD
> Alameda, 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.
>

[[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] pattern matching

2013-01-07 Thread Data Analytics Corp.

Hi,

I have a simple question.  Suppose I have a string "x$Expensive". I want 
to find the position of the $ in this string; i.e., I want a function 
that returns 2.  I tried grep, regexpr, etc with no luck, unless I'm 
just using them incorrectly.  Any suggestions?


Thanks,

Walt



Walter R. Paczkowski, Ph.D.
Data Analytics Corp.
44 Hamilton Lane
Plainsboro, NJ 08536

(V) 609-936-8999
(F) 609-936-3733
w...@dataanalyticscorp.com
www.dataanalyticscorp.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] try()-function does not catch error in BATCH-job if Matrix is loaded

2013-01-07 Thread luke-tierney

This is due to long-staning issue in methods internals, which are
involved because loading Matrix shadows base::mean with Matrix::mean.
A work-around has been in place in R_devel for some time; a proper fix
may come at some point in the future. So if your real code doesn't
need the moficied mean from Matrix you can use base::mean; otherwise
you will have to use an R-devel snapshot.

Best,

luke

On Mon, 7 Jan 2013, Sarah Brockhaus wrote:


Hello,

In my simulation I use the try()-function to catch possible errors when 
fitting models. I run the simulationon a Linux-server using the command " R 
CMD BATCH nameOfFile.R &".  When executing the code as batch-job I get the 
problem that the execution is halted without giving an error message. But 
when I run the code interactivly the error is catched by try() as it would be 
expected.


The problem is somewhat strange as it only occurs when the code is executed 
as a batch-job and when the package "Matrix" is loaded.


I wrote a small example reproducing the error. (In my code the error occurs 
in mgcv:::fixDependence, which looks like the code I'm using below in order 
to get a small reproducible example. I realized that the code  makes no 
sense...)


##
library(Matrix)

R <- matrix(abs(rnorm(25)), 5, 5)
r0 <- r <- nrow(R)

# while-loop produces error that should be catched by the function try()
try(
while (mean(R[r0:r, r0:r]) > 0) r0 <- r0 - 1
)

# so "Hello" should be printed
print("Hello")
##

I use R version 2.15.1 (2012-06-22) on a x86_64-pc-linux-gnu (64-bit) 
platform.


I would be grateful for help.

Best regards,
Sarah Brockhaus

PHD-student
Department of Statistics
University of Munich

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



--
Luke Tierney
Chair, Statistics and Actuarial Science
Ralph E. Wareham Professor of Mathematical Sciences
University of Iowa  Phone: 319-335-3386
Department of Statistics andFax:   319-335-3017
   Actuarial Science
241 Schaeffer Hall  email:   luke-tier...@uiowa.edu
Iowa City, IA 52242 WWW:  http://www.stat.uiowa.edu

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] random effects model

2013-01-07 Thread arun
HI,

Regarding question:2) Have you checked summary(m1)?
data(seizure)
 ## Diggle, Liang, and Zeger (1994) pp166-168, compare Table 8.10
 seiz.l <- reshape(seizure,
   varying=list(c("base","y1", "y2", "y3", "y4")),
   v.names="y", times=0:4, direction="long")
 seiz.l <- seiz.l[order(seiz.l$id, seiz.l$time),]
 seiz.l$t <- ifelse(seiz.l$time == 0, 8, 2)
 seiz.l$x <- ifelse(seiz.l$time == 0, 0, 1)
 m1 <- geese(y ~ offset(log(t)) + x + trt + x:trt, id = id,
 data=seiz.l, corstr="exch", family=poisson)
 summary(m1)
#Call:
#geese(formula = y ~ offset(log(t)) + x + trt + x:trt, id = id, 
 #   data = seiz.l, family = poisson, corstr = "exch")
#
#Mean Model:
# Mean Link: log 
# Variance to Mean Relation: poisson 

 #Coefficients:
  #     estimate    san.se   wald p   
#(Intercept)  1.34760922 0.1573571 73.3423807 0.000
#x    0.11183602 0.1159304  0.9306116 0.3347040
#trt  0.02753449 0.2217878  0.0154127 0.9011982
#x:trt   -0.10472579 0.2134448  0.2407334 0.6236769

#p is the colname with the p values for the Coefficients
summary(m1)$mean["p"]
#    p
#(Intercept) 0.000
#x   0.3347040
#trt 0.9011982
#x:trt   0.6236769


You didn't say specifically whether you got NA's in example data or your actual 
data.  I am getting the p-values in R 2.15. 


1) I think it should work with binary response variable.
Another example in the same documentation:
 data(ohio)
str(ohio)
#'data.frame':    2148 obs. of  4 variables:
# $ resp : int  0 0 0 0 0 0 0 0 0 0 ...
# $ id   : int  0 0 0 0 1 1 1 1 2 2 ...
# $ age  : int  -2 -1 0 1 -2 -1 0 1 -2 -1 ...
# $ smoke: int  0 0 0 0 0 0 0 0 0 0 ...
It is not even factors


 fit <- geese(resp ~ age + smoke + age:smoke, id=id, data=ohio,
  family=binomial, corstr="exch", scale.fix=TRUE)
 summary(fit)$mean["p"]
 #    p
#(Intercept) 0.
#age 0.01523698
#smoke   0.09478252
#age:smoke   0.42234200


# also tested after converting to factors

ohio1<-ohio
ohio1$smoke<-factor(ohio1$smoke)
 ohio1$age<-factor(ohio1$age)


str(ohio1)
#'data.frame':    2148 obs. of  4 variables:
# $ resp : int  0 0 0 0 0 0 0 0 0 0 ...
# $ id   : int  0 0 0 0 1 1 1 1 2 2 ...
# $ age  : Factor w/ 4 levels "-2","-1","0",..: 1 2 3 4 1 2 3 4 1 2 ...
# $ smoke: Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
fit1<-geese(resp~age+smoke+age:smoke,id=id,data=ohio1,family=binomial,corstr="exch",scale.fix=TRUE)
 summary(fit1)$mean["p"]
#  p
#(Intercept)  0.
#age-1    0.60555454
#age0 0.45322698
#age1 0.01187725
#smoke1   0.86262269
#age-1:smoke1 0.17239050
#age0:smoke1  0.32223942
#age1:smoke1  0.36686706


A.K.

- Original Message -
From: rex2013 
To: r-help@r-project.org
Cc: 
Sent: Monday, January 7, 2013 6:15 AM
Subject: Re: [R] random effects model

Hi A.K

Below is the comment I get, not sure why.

BP.sub3 is the stacked data without the missing values.

BP.geese3 <- geese(HiBP~time*MaternalAge,data=BP.sub3,id=CODEA,
family=binomial, corstr="unstructured", na.action=na.omit)Error in
`contrasts<-`(`*tmp*`, value = contr.funs[1 + isOF[nn]]) :
  contrasts can be applied only to factors with 2 or more levels

Even though age has 3 levels; time has 14 years & 21 years; HIBP is a
binary response outcome.

2) When you mentioned summary(m1)$mean["p"] what did the p mean? i
used this in one of the gee command, it produced NA as answer?

Many thanks



On Mon, Jan 7, 2013 at 5:26 AM, arun kirshna [via R] <
ml-node+s789695n4654795...@n4.nabble.com> wrote:

> Hi,
>
> I am  not very familiar with the geese/geeglm().  Is it from
> library(geepack)?
> Regarding your question:
> "
> Can you tell me if I can use the geese or geeglm function with this data
> eg: : HIBP~ time* Age
> Here age is a factor with 3 levels, time: 2 levels, HIBP = yes/no.
>
> From your original data:
> BP_2b<-read.csv("BP_2b.csv",sep="\t")
> head(BP_2b,2)
> #  CODEA Sex MaternalAge Education Birthplace AggScore IntScore Obese14
> #1     1  NA           3         4          1       NA       NA      NA
> #2     3   2           3         3          1        0        0       0
>  # Overweight14 Overweight21 Obese21 hibp14 hibp21
> #1           NA           NA      NA     NA     NA
> #2            0            1       0      0      0
>
> If I understand your new classification:
> BP.stacknormal<- subset(BP_2b,Obese14==0 & Overweight14==0 & Obese21==0 &
> Overweight21==0)
> BP.stackObese <- subset(BP_2b,(Obese14==1& Overweight14==0 &
> Obese14==1&Overweight14==1)|(Obese14==1&Overweight14==1 & Obese21==1 &
> Overweight21==0)|(Obese14==1&Overweight14==0 & Obese21==0 &
> Overweight21==0)|(Obese14==0 & Overweight14==0 & Obese21==1 &
> Overweight21==0)|(Obese14==0 & Overweight14==0 & Obese21==1 &
> Ov

Re: [R] random effects model

2013-01-07 Thread arun
HI,

Just to add:
fit3<-geese(hibp~MaternalAge*time,id=CODEA,data=BP.stack5,family=binomial,corstr="exch",scale.fix=TRUE)
 #works
 summary(fit3)$mean["p"]
# p
#(Intercept) 0.
#MaternalAge4    0.49099242
#MaternalAge5    0.04686295
#time21  0.86164351
#MaternalAge4:time21 0.59258221
#MaternalAge5:time21 0.79909832

fit4<-geese(hibp~MaternalAge*time,id=CODEA,data=BP.stack5,family=binomial,corstr="unstructured",scale.fix=TRUE)
  #when the correlation structure is changed to "unstructured"
#Error in `contrasts<-`(`*tmp*`, value = contr.funs[1 + isOF[nn]]) : 
 # contrasts can be applied only to factors with 2 or more levels
#In addition: Warning message:
#In is.na(rows) : is.na() applied to non-(list or vector) of type 'NULL'


Though, it works with data(Ohio)

fit1<-geese(resp~age+smoke+age:smoke,id=id,data=ohio1,family=binomial,corstr="unstructured",scale.fix=TRUE)
 summary(fit1)$mean["p"]
#  p
#(Intercept)  0.
#age-1    0.60555454
#age0 0.45322698
#age1 0.01187725
#smoke1   0.86262269
#age-1:smoke1 0.17239050
#age0:smoke1  0.32223942
#age1:smoke1  0.36686706



By checking:
 with(BP.stack5,table(MaternalAge,time))
#   time
#MaternalAge   14   21
  #    3 1104  864
   #   4  875  667
    #     5   67   53 #less number of observations


 BP.stack6 <- BP.stack5[order(BP.stack5$CODEA, BP.stack5$time),]
 head(BP.stack6)  # very few IDs with  MaternalAge==5
#   X CODEA Sex MaternalAge Education Birthplace AggScore IntScore
#1493 3.1 3   2   3 3  1    0    0
#3202 3.2 3   2   3 3  1    0    0
#1306 7.1 7   2   4 6  1    0    0
#3064 7.2 7   2   4 6  1    0    0
#1    8.1 8   2   4 4  1    0    0
#2047 8.2 8   2   4 4  1    0    0
 #     Categ time Obese Overweight hibp
#1493 Overweight   14 0  0    0
#3202 Overweight   21 0  1    0
#1306  Obese   14 0  0    0
#3064  Obese   21 1  1    0
#1    Normal   14 0  0    0
#2047 Normal   21 0  0    0
BP.stack7<-BP.stack6[BP.stack6$MaternalAge!=5,]
 BP.stack7$MaternalAge<-factor(as.numeric(as.character(BP.stack7$MaternalAge)
 
fit5<-geese(hibp~MaternalAge*time,id=CODEA,data=BP.stack7,family=binomial,corstr="unstructured",scale.fix=TRUE)
#Error in `contrasts<-`(`*tmp*`, value = contr.funs[1 + isOF[nn]]) : 
 # contrasts can be applied only to factors with 2 or more levels

 with(BP.stack7,table(MaternalAge,time))  #It looks like the combinations are 
still there
#   time
#MaternalAge   14   21
 #     3 1104  864
   #   4  875  667

It works also with corstr="ar1".   Why do you gave the option "unstructured"? 
A.K.






- Original Message -
From: rex2013 
To: r-help@r-project.org
Cc: 
Sent: Monday, January 7, 2013 6:15 AM
Subject: Re: [R] random effects model

Hi A.K

Below is the comment I get, not sure why.

BP.sub3 is the stacked data without the missing values.

BP.geese3 <- geese(HiBP~time*MaternalAge,data=BP.sub3,id=CODEA,
family=binomial, corstr="unstructured", na.action=na.omit)Error in
`contrasts<-`(`*tmp*`, value = contr.funs[1 + isOF[nn]]) :
  contrasts can be applied only to factors with 2 or more levels

Even though age has 3 levels; time has 14 years & 21 years; HIBP is a
binary response outcome.

2) When you mentioned summary(m1)$mean["p"] what did the p mean? i
used this in one of the gee command, it produced NA as answer?

Many thanks



On Mon, Jan 7, 2013 at 5:26 AM, arun kirshna [via R] <
ml-node+s789695n4654795...@n4.nabble.com> wrote:

> Hi,
>
> I am  not very familiar with the geese/geeglm().  Is it from
> library(geepack)?
> Regarding your question:
> "
> Can you tell me if I can use the geese or geeglm function with this data
> eg: : HIBP~ time* Age
> Here age is a factor with 3 levels, time: 2 levels, HIBP = yes/no.
>
> From your original data:
> BP_2b<-read.csv("BP_2b.csv",sep="\t")
> head(BP_2b,2)
> #  CODEA Sex MaternalAge Education Birthplace AggScore IntScore Obese14
> #1     1  NA           3         4          1       NA       NA      NA
> #2     3   2           3         3          1        0        0       0
>  # Overweight14 Overweight21 Obese21 hibp14 hibp21
> #1           NA           NA      NA     NA     NA
> #2            0            1       0      0      0
>
> If I understand your new classification:
> BP.stacknormal<- subset(BP_2b,Obese14==0 & Overweight14==0 & Obese21==0 &
> Overweight21==0)
> BP.stackObese <- subset(BP_2b,(Obese14==1& Overweight14==0 &
> Obese14==1&Overweight14==1)|(Obese14==1&Overweight14==1 & Obese21==1 &
> Overweight21==0)|(Obese14==1&Overweight14==0 & Obese21==0 &
> Overweight21==0)|(Obese14==0 & Overweight14==0 & Obese21==1 &
> Overweight21==0)|(Obese

Re: [R] Changing mtext direction, or using text for the margin?

2013-01-07 Thread Michael Rennie


HA! Works perfectly!!! Thanks so much ilai.

ilai wrote:


On Mon, Jan 7, 2013 at 11:52 AM, David Winsemius 
mailto:dwinsem...@comcast.net>> wrote:



If you are going that route you may want to look at the gridBase
package.

Yes for mixing base and grid graphics but IMHO overkill here. 
Replacing the last mtext line with


grid::grid.text('dependent B', 0.985 , 0.5 , rot = 270)

should take care of it (leaving exact placement and justification to 
the OP)


cheers

-- 
David.


Mike

David Winsemius wrote:


On Jan 7, 2013, at 8:06 AM, Michael Rennie wrote:


Any thoughts on what that dirty hack might be or any
leads on where to start? Perhaps a whole new plot
region in the margin or something? Is that even
possible? I'm having a difficult time imagining how I
can do this.


The `text` function accepts the 'xpd' argument for `par`.

?par
?text

(Similar but not exactly the same as R-FAQ "7.27 How can I
create rotated axis labels?")

I'm not sure what sort of "dirt" Uwe was expecting, but
I'm guessing he meant that you probably were not going to
be getting exactly what you wanted with the first attempt
and were going to need to adjust the locations "by hand"
after digging around in the par documentation.



David Winsemius, MD
Alameda, 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.




__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] pattern matching

2013-01-07 Thread Marc Schwartz
On Jan 7, 2013, at 3:22 PM, Data Analytics Corp.  
wrote:

> Hi,
> 
> I have a simple question.  Suppose I have a string "x$Expensive". I want to 
> find the position of the $ in this string; i.e., I want a function that 
> returns 2.  I tried grep, regexpr, etc with no luck, unless I'm just using 
> them incorrectly.  Any suggestions?
> 
> Thanks,
> 
> Walt



The problem with this specific example is that '$' is a metacharacter in 
regular expressions, so you have to escape it. For example:

> regexpr("\\$", "x$Expensive")
[1] 2
attr(,"match.length")
[1] 1
attr(,"useBytes")
[1] TRUE



See ?regex for more information and if appropriate, consider gregexpr():

> gregexpr("\\$", "x$Expensive$MoreText")
[[1]]
[1]  2 12
attr(,"match.length")
[1] 1 1
attr(,"useBytes")
[1] TRUE


Regards,

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] pattern matching

2013-01-07 Thread William Dunlap
"$" has a special meaning (end-of-string) in regular expressions, so you can 
either escape it with "\\" or not use regular expressions in regexpr():

> regexpr("\\$", "x$Expensive")
[1] 2
attr(,"match.length")
[1] 1
attr(,"useBytes")
[1] TRUE
> regexpr("$", "x$Expensive", fixed=TRUE)
[1] 2
attr(,"match.length")
[1] 1
attr(,"useBytes")
[1] TRUE


Bill Dunlap
Spotfire, TIBCO Software
wdunlap tibco.com


> -Original Message-
> From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On 
> Behalf
> Of Data Analytics Corp.
> Sent: Monday, January 07, 2013 1:22 PM
> To: r-help@R-project.org
> Subject: [R] pattern matching
> 
> Hi,
> 
> I have a simple question.  Suppose I have a string "x$Expensive". I want
> to find the position of the $ in this string; i.e., I want a function
> that returns 2.  I tried grep, regexpr, etc with no luck, unless I'm
> just using them incorrectly.  Any suggestions?
> 
> Thanks,
> 
> Walt
> 
> 
> 
> Walter R. Paczkowski, Ph.D.
> Data Analytics Corp.
> 44 Hamilton Lane
> Plainsboro, NJ 08536
> 
> (V) 609-936-8999
> (F) 609-936-3733
> w...@dataanalyticscorp.com
> www.dataanalyticscorp.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] pattern matching

2013-01-07 Thread arun
HI,
str1<-"x$Expensive"
regexpr("\\$",str1)[1]
#[1] 2
 str2<-"x$Exp$Expression"
unlist(gregexpr("\\$",str2))
#[1] 2 6


A.K.




- Original Message -
From: Data Analytics Corp. 
To: "r-help@R-project.org" 
Cc: 
Sent: Monday, January 7, 2013 4:22 PM
Subject: [R] pattern matching

Hi,

I have a simple question.  Suppose I have a string "x$Expensive". I want to 
find the position of the $ in this string; i.e., I want a function that returns 
2.  I tried grep, regexpr, etc with no luck, unless I'm just using them 
incorrectly.  Any suggestions?

Thanks,

Walt



Walter R. Paczkowski, Ph.D.
Data Analytics Corp.
44 Hamilton Lane
Plainsboro, NJ 08536

(V) 609-936-8999
(F) 609-936-3733
w...@dataanalyticscorp.com
www.dataanalyticscorp.com

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


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


[R] Amelia algorithm

2013-01-07 Thread Martin
Dear all.

First of all, my english isn't verry good, but I hope I can convey my concern.
I've a general question about the Amelia algorithm. I'm no mathematician or
statistician, but I had to use R and impute and analyse some data, and Amelia
showed results that fitted my expectations. I'll have to defend my choice soon,
but I haven't totally grasped what Amelia does. I'm particularly interested in a
simple as possible explanation in how Amelia imputation works. I've read that it
uses a bootstrapping-based algorithm, but how does it chose the values?
The data had mainly value >0 (chemical concentrations, water temperature and
pH-value).

Regards Martin

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] checker/chequer board pattern as a colour

2013-01-07 Thread Greg Snow
Not as a simple color, but you can use the rasterImage function to add a
rectangle with a checkerboard pattern to a an existing plot, the examples
for the function use a simple checkerboard pattern.

Be careful of overusing something like this, too many patterned areas in a
plot can be more distracting and distorting than useful.


On Mon, Jan 7, 2013 at 8:00 AM, Federico Calboli
wrote:

> Hi All,
>
> is there a reasonably simple way of using a black and white
> chequer/checker board pattern as a colour:
>
> barplot(mydata, col = c('red', 'blue' 'checkerboard'))
>
> ?
>
> BW
>
> F
>
> --
> Federico C. F. Calboli
> Neuroepidemiology and Ageing Research
> Imperial College, St. Mary's Campus
> Norfolk Place, London W2 1PG
>
> Tel +44 (0)20 75941602   Fax +44 (0)20 75943193
>
> f.calboli [.a.t] imperial.ac.uk
> f.calboli [.a.t] 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.
>



-- 
Gregory (Greg) L. Snow Ph.D.
538...@gmail.com

[[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] Changing mtext direction, or using text for the margin?

2013-01-07 Thread Greg Snow
And look at the grconvertX and grconvertY functions for possibly ways to
compute positioning information when using text with base graphics.


On Mon, Jan 7, 2013 at 11:35 AM, David Winsemius wrote:

>
> On Jan 7, 2013, at 8:06 AM, Michael Rennie wrote:
>
>
>> Any thoughts on what that dirty hack might be or any leads on where to
>> start? Perhaps a whole new plot region in the margin or something? Is that
>> even possible? I'm having a difficult time imagining how I can do this.
>>
>
> The `text` function accepts the 'xpd' argument for `par`.
>
> ?par
> ?text
>
> (Similar but not exactly the same as R-FAQ "7.27 How can I create rotated
> axis labels?")
>
> I'm not sure what sort of "dirt" Uwe was expecting, but I'm guessing he
> meant that you probably were not going to be getting exactly what you
> wanted with the first attempt and were going to need to adjust the
> locations "by hand" after digging around in the par documentation.
>
> --
> David.
>
>  Mike
>>
>> Uwe Ligges wrote:
>>
>>>
>>>
>>> On 07.01.2013 07:00, Michael Rennie wrote:
>>>

 Hi all,

 I have read through the archives, but can't find a solution to this
 problem.

 I need the text direction on "dependent B", plotted in margin 4, to go
 top to bottom (opposite what it is now). Here's some sample code:

 #plot with mtext example

 par(mgp = c(2,1,0), mfrow=c(2,2), las=1, mar=c(2,2,2,2), omi=
 c(0.5,0.2,0,0.2))

 a<-1:10
 b<-7:16
 c<-21:30

 plot(a~b, ylab="", xlab="", xaxt="n")
 axis(1, at=b, labels=FALSE, tick=TRUE)
 plot(c~b, ylab="", xlab="", xaxt="n", yaxt="n")
 axis(1, at=b, labels=FALSE, tick=TRUE)
 axis(4, at=c(22,24,26,28,30), labels=TRUE, tick=TRUE)
 plot(a~b, ylab="")
 plot(c~b, ylab="", yaxt="n")
 axis(4, at=c(22,24,26,28,30), labels=TRUE, tick=TRUE)

 mtext(side=1, "independent", outer=TRUE, line=1, padj=1)
 mtext(side=2, "dependent A", las=0, outer=TRUE, line=0.25)
 mtext(side=4, "dependent B", las=0, outer=TRUE, line=0.25)

 I have seen an example in help where i can use text and (srt) to
 manipulate this with a single panel plot, but not for a multi-panel
 example. Any help would be greatly appreciated.

>>>
>>>
>>> You cannot do that with mtext, you rather need a dirty hack with text(),
>>> I believe.
>>>
>>> Best,
>>> Uwe Ligges
>>>
>>>
>>>
>>>  Cheers,

 Mike

 __**
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/**listinfo/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.
>>
>
> David Winsemius, MD
> Alameda, 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.
>



-- 
Gregory (Greg) L. Snow Ph.D.
538...@gmail.com

[[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] checker/chequer board pattern as a colour

2013-01-07 Thread Federico Calboli
On 7 Jan 2013, at 21:53, Greg Snow <538...@gmail.com> wrote:

> Not as a simple color, but you can use the rasterImage function to add a 
> rectangle with a checkerboard pattern to a an existing plot, the examples for 
> the function use a simple checkerboard pattern.  
> 
> Be careful of overusing something like this, too many patterned areas in a 
> plot can be more distracting and distorting than useful.

thanks, but that sounds way too much hassle than it's worth

BW


F




> 
> 
> On Mon, Jan 7, 2013 at 8:00 AM, Federico Calboli  
> wrote:
> Hi All,
> 
> is there a reasonably simple way of using a black and white chequer/checker 
> board pattern as a colour:
> 
> barplot(mydata, col = c('red', 'blue' 'checkerboard'))
> 
> ?
> 
> BW
> 
> F
> 
> --
> Federico C. F. Calboli
> Neuroepidemiology and Ageing Research
> Imperial College, St. Mary's Campus
> Norfolk Place, London W2 1PG
> 
> Tel +44 (0)20 75941602   Fax +44 (0)20 75943193
> 
> f.calboli [.a.t] imperial.ac.uk
> f.calboli [.a.t] 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.
> 
> 
> 
> -- 
> Gregory (Greg) L. Snow Ph.D.
> 538...@gmail.com


--
Federico C. F. Calboli
Neuroepidemiology and Ageing Research
Imperial College, St. Mary's Campus
Norfolk Place, London W2 1PG

Tel +44 (0)20 75941602   Fax +44 (0)20 75943193

f.calboli [.a.t] imperial.ac.uk
f.calboli [.a.t] 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] Changing mtext direction, or using text for the margin?

2013-01-07 Thread David L Carlson
Or use the "ugly hack" with a bit of locator() trial and error.
Replace the last mtext() call with

text(17.9, 31.9, "dependent B", srt=-90, xpd=NA, cex=1.25)

--
David L Carlson
Associate Professor of Anthropology
Texas A&M University
College Station, TX 77843-4352


> -Original Message-
> From: r-help-boun...@r-project.org [mailto:r-help-bounces@r-
> project.org] On Behalf Of Greg Snow
> Sent: Monday, January 07, 2013 3:59 PM
> To: David Winsemius
> Cc: r-help; Uwe Ligges; Michael Rennie
> Subject: Re: [R] Changing mtext direction, or using text for the
> margin?
> 
> And look at the grconvertX and grconvertY functions for possibly ways
> to
> compute positioning information when using text with base graphics.
> 
> 
> On Mon, Jan 7, 2013 at 11:35 AM, David Winsemius
> wrote:
> 
> >
> > On Jan 7, 2013, at 8:06 AM, Michael Rennie wrote:
> >
> >
> >> Any thoughts on what that dirty hack might be or any leads on where
> to
> >> start? Perhaps a whole new plot region in the margin or something?
> Is that
> >> even possible? I'm having a difficult time imagining how I can do
> this.
> >>
> >
> > The `text` function accepts the 'xpd' argument for `par`.
> >
> > ?par
> > ?text
> >
> > (Similar but not exactly the same as R-FAQ "7.27 How can I create
> rotated
> > axis labels?")
> >
> > I'm not sure what sort of "dirt" Uwe was expecting, but I'm guessing
> he
> > meant that you probably were not going to be getting exactly what you
> > wanted with the first attempt and were going to need to adjust the
> > locations "by hand" after digging around in the par documentation.
> >
> > --
> > David.
> >
> >  Mike
> >>
> >> Uwe Ligges wrote:
> >>
> >>>
> >>>
> >>> On 07.01.2013 07:00, Michael Rennie wrote:
> >>>
> 
>  Hi all,
> 
>  I have read through the archives, but can't find a solution to
> this
>  problem.
> 
>  I need the text direction on "dependent B", plotted in margin 4,
> to go
>  top to bottom (opposite what it is now). Here's some sample code:
> 
>  #plot with mtext example
> 
>  par(mgp = c(2,1,0), mfrow=c(2,2), las=1, mar=c(2,2,2,2), omi=
>  c(0.5,0.2,0,0.2))
> 
>  a<-1:10
>  b<-7:16
>  c<-21:30
> 
>  plot(a~b, ylab="", xlab="", xaxt="n")
>  axis(1, at=b, labels=FALSE, tick=TRUE)
>  plot(c~b, ylab="", xlab="", xaxt="n", yaxt="n")
>  axis(1, at=b, labels=FALSE, tick=TRUE)
>  axis(4, at=c(22,24,26,28,30), labels=TRUE, tick=TRUE)
>  plot(a~b, ylab="")
>  plot(c~b, ylab="", yaxt="n")
>  axis(4, at=c(22,24,26,28,30), labels=TRUE, tick=TRUE)
> 
>  mtext(side=1, "independent", outer=TRUE, line=1, padj=1)
>  mtext(side=2, "dependent A", las=0, outer=TRUE, line=0.25)
>  mtext(side=4, "dependent B", las=0, outer=TRUE, line=0.25)
> 
>  I have seen an example in help where i can use text and (srt) to
>  manipulate this with a single panel plot, but not for a multi-
> panel
>  example. Any help would be greatly appreciated.
> 
> >>>
> >>>
> >>> You cannot do that with mtext, you rather need a dirty hack with
> text(),
> >>> I believe.
> >>>
> >>> Best,
> >>> Uwe Ligges
> >>>
> >>>
> >>>
> >>>  Cheers,
> 
>  Mike
> 
>  __**
>  R-help@r-project.org mailing list
>  https://stat.ethz.ch/mailman/**listinfo/r-
> help
>  PLEASE do read the posting guide
>  http://www.R-project.org/**posting-guide.html 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.
> >>
> >
> > David Winsemius, MD
> > Alameda, 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.
> >
> 
> 
> 
> --
> Gregory (Greg) L. Snow Ph.D.
> 538...@gmail.com
> 
>   [[alternative HTML version deleted]]
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-
> guide.html
> and provide commented, minimal, self-contained, reproducible code.

___

Re: [R] how to aggregate T-test result in an elegant way?

2013-01-07 Thread Yao He
Yes, thanks a lot for your help!

Regards

2013/1/8 arun :
> Hi Yao,
>
> You could also have the results in a wide format:
> res<-do.call(rbind,lapply(lapply(split(b,b$variable),function(x) 
> t.test(x$value[x$O2=="13%"],x$value[x$O2=="21%"])),function(x) 
> data.frame(mean13=x$estimate[1],mean21=x$estimate[2],p.value=x$p.value,CILow=x$conf.int[1],CIHigh=x$conf.int[2])))
>  res
> #  mean13   mean21p.value CILowCIHigh
> #EMW 14.35000 17.68000 0.09355374 -7.682686  1.022686
> #EW.17.5 42.87000 45.97333 0.17464018 -9.265622  3.058955
> #EW.INCU 49.61333 47.08333 0.43689727 -7.119234 12.179234
> A.K.
>
>
>
>
> - Original Message -
> From: Yao He 
> To: arun 
> Cc: R help 
> Sent: Monday, January 7, 2013 10:57 AM
> Subject: Re: [R] how to aggregate T-test result in an elegant way?
>
> Hi,arun
>
> Yes , I just want to do the t.test
> I think maybe  it is not necessary to generate a 3D array from the raw
> data.frame by acast() at first
>
> Thanks a lot
>
> 2013/1/7 arun :
>> Hi Yao,
>>
>> It's okay.
>>
>> How did you generate the 3 D array?
>> Using ?acast()
>>
>> I am not sure I understand your question "
>>
>> if you meet a t-test task as I described  , is that generate a
>> high-dimension array  a good way ?"
>>
>> Do you want to do the t-test in the melt dataset?
>>
>> b<- read.table(text="
>> IDO2variablevalue
>> 1TWF2H513% EW.INCU49.38
>> 2TWF2H613% EW.INCU48.02
>> 3TWF2H1913%EW.INCU51.44
>> 280TWF2H10113% EW.17.542.26
>> 281TWF2H10513%EW.17.543.52
>> 282TWF2H10613% EW.17.542.83
>> 472TWF2N10221% EW.17.545.97
>> 473TWF2N10421%EW.17.5 43.32
>> 474TWF2N10621% EW.17.548.63
>> 689TWF2N221% EMW19.57
>> 690TWF2N621%EMW18.07
>> 691TWF2N1021%EMW15.4
>> 491TWF2H513%EMW15.61
>> 492TWF2H613%EMW13.41
>> 493TWF2H1913%EMW14.03
>> 199TWF2N221%EW.INCU48.69
>> 200TWF2N621%EW.INCU50.52
>> 201TWF2N1021%EW.INCU42.04
>> ",sep="",header=TRUE,stringsAsFactors=FALSE)
>>  res<-lapply(lapply(split(b,b$variable),function(x) 
>> t.test(x$value[x$O2=="13%"],x$value[x$O2=="21%"])),function(x) 
>> data.frame(mean=x$estimate,p.value=x$p.value))
>> res1<-do.call(rbind,res)
>> row.names(res1)[grep("mean of 
>> x",row.names(res1))]<-gsub("(.*\\.).*$","\\113%",row.names(res1)[grep("mean 
>> of x",row.names(res1))])
>> row.names(res1)[grep("mean of 
>> y",row.names(res1))]<-gsub("(.*\\.).*$","\\121%",row.names(res1)[grep("mean 
>> of y",row.names(res1))])
>> res1
>> #meanp.value
>> #EMW.13% 14.35000 0.09355374
>> #EMW.21% 17.68000 0.09355374
>> #EW.17.5.13% 42.87000 0.17464018
>> #EW.17.5.21% 45.97333 0.17464018
>> #EW.INCU.13% 49.61333 0.43689727
>> #EW.INCU.21% 47.08333 0.43689727
>>
>> A.K.
>>
>>
>>
>> - Original Message -
>> From: Yao He 
>> To: arun 
>> Cc: R help 
>> Sent: Monday, January 7, 2013 4:00 AM
>> Subject: Re: [R] how to aggregate T-test result in an elegant way?
>>
>> Hi, arun
>> I'm so sorry for that isn't helpful.
>> One of question is that I don't know how  to subset a small part as it
>> is a 3-dimension array so I just show the structure of that.
>> I tried  dput()  to a file , then what should I do for subsetting it?
>>
>> Another question is :
>> My rawdata is a "melt" dataframe like that:
>> IIDO2variablevalue
>> 1TWF2H513% EW.INCU49.38
>> 2TWF2H613% EW.INCU48.02
>> 3TWF2H1913% EW.INCU51.44
>> 280TWF2H10113% EW.17.542.26
>> 281TWF2H10513% EW.17.5 43.52
>> 282TWF2H10613% EW.17.542.83
>> 472TWF2N10221% EW.17.545.97
>> 473TWF2N10421% EW.17.5 43.32
>> 474TWF2N10621% EW.17.548.63
>> 689TWF2N221%  EMW19.57
>> 690TWF2N621% EMW18.07
>> 691TWF2N1021%EMW15.4
>> 491TWF2H5 13%EMW15.61
>> 492TWF2H613%EMW13.41
>> 493TWF2H1913%EMW14.03
>> 199TWF2N221%EW.INCU48.69
>> 200TWF2N621%EW.INCU50.52
>> 201TWF2N1021%EW.INCU42.04
>>
>> if you meet a t-test task as I described  , is that generate a
>> high-dimension array  a good way ?
>> Thank you!
>>
>> Yao He
>> 2013/1/7 arun :
>>> HI,
>>> I tried to create an example dataset (as you didn't provide the data).
>>> set.seed(25)
>>> a<-array(sample(1:50,60,replace=TRUE),dim=c(2,10,3))
>>> dimnames(a)[[1]]<-c("13%","21%")
>>> dimnames(a)[[2]]<-paste("TWF2H",101:110,sep="")
>>> dimnames(a)[[3]]<-c("EW.INCU","EW.17.5","EMW")
>>>
>>>
>>> str(a)
>>> # int [1:2, 1:10, 1:3] 21 35 8 45 7 50 32 17 4 15 ...
>>>  #- attr(*, "dimnames")=List of 3
>>>   #..$ : chr [1:2] "13%" "21%"
>>>   .#.$ : chr [1:10] "TWF2H101" "TWF2H102" "TWF2H103

Re: [R] Amelia algorithm

2013-01-07 Thread zGreenfelder
On Mon, Jan 7, 2013 at 4:29 PM, Martin  wrote:
> Dear all.
>
> First of all, my english isn't verry good, but I hope I can convey my concern.
> I've a general question about the Amelia algorithm. I'm no mathematician or
> statistician, but I had to use R and impute and analyse some data, and Amelia
> showed results that fitted my expectations. I'll have to defend my choice 
> soon,
> but I haven't totally grasped what Amelia does. I'm particularly interested 
> in a
> simple as possible explanation in how Amelia imputation works. I've read that 
> it
> uses a bootstrapping-based algorithm, but how does it chose the values?
> The data had mainly value >0 (chemical concentrations, water temperature and
> pH-value).
>
> Regards Martin

I'm pretty new here, but a quick google search suggests that perhaps
http://gking.harvard.edu/amelia (and maybe google translate)
might have some decent pointers for you.

I poked at the documentation from that site (a pdf file), and it's
quite intense on
the mathematics, you may get more from it than I could.  there's also
a link there
to a separate, for thisspecific package/algorithm  (it seems to be
called Amel ia II,
 I'm assuming this is the same that you used, if I'm off .. sorry about that)

HTH
-- 
Even the Magic 8 ball has an opinion on email clients: Outlook not so good.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] random effects model

2013-01-07 Thread arun
HI,
BP.stack5 is the one without missing values.
na.omit().  Otherwise, I have to use the option na.action=.. in the 
?geese() statement

You need to read about the correlation structures.  IN unstructured option, 
more number of parameters needs to be estimated,  In repeated measures design, 
when the underlying structure is not known, it would be better to compare using 
different options (exchangeable is similar to compound symmetry) and select the 
one which provide the least value for AIC or BIC.  Have a look at 

http://stats.stackexchange.com/questions/21771/how-to-perform-model-selection-in-gee-in-r
It's not clear to me  "reference to write about missing values".    
A.K.




- Original Message -
From: Usha Gurunathan 
To: arun 
Cc: 
Sent: Monday, January 7, 2013 6:12 PM
Subject: Re: [R] random effects model

Hi AK

2)I shall try putting exch. and check when I get home. Btw, what is
BP.stack5? is it with missing values or only complete cases?

I guess I am still not clear about the unstructured and exchangeable
options, as in which one is better.

1)Rgding the summary(p): NA thing, I tried putting one of my gee equation.

Can you suggest me a reference to write about" missing values and the
implications for my results"

Thanks.



On 1/8/13, arun  wrote:
> HI,
>
> Just to add:
> fit3<-geese(hibp~MaternalAge*time,id=CODEA,data=BP.stack5,family=binomial,corstr="exch",scale.fix=TRUE)
> #works
>  summary(fit3)$mean["p"]
> #                             p
> #(Intercept)         0.
> #MaternalAge4        0.49099242
> #MaternalAge5        0.04686295
> #time21              0.86164351
> #MaternalAge4:time21 0.59258221
> #MaternalAge5:time21 0.79909832
>
> fit4<-geese(hibp~MaternalAge*time,id=CODEA,data=BP.stack5,family=binomial,corstr="unstructured",scale.fix=TRUE)
> #when the correlation structure is changed to "unstructured"
> #Error in `contrasts<-`(`*tmp*`, value = contr.funs[1 + isOF[nn]]) :
>  # contrasts can be applied only to factors with 2 or more levels
> #In addition: Warning message:
> #In is.na(rows) : is.na() applied to non-(list or vector) of type 'NULL'
>
>
> Though, it works with data(Ohio)
>
> fit1<-geese(resp~age+smoke+age:smoke,id=id,data=ohio1,family=binomial,corstr="unstructured",scale.fix=TRUE)
>  summary(fit1)$mean["p"]
> #                      p
> #(Intercept)  0.
> #age-1        0.60555454
> #age0         0.45322698
> #age1         0.01187725
> #smoke1       0.86262269
> #age-1:smoke1 0.17239050
> #age0:smoke1  0.32223942
> #age1:smoke1  0.36686706
>
>
>
> By checking:
>  with(BP.stack5,table(MaternalAge,time))
> #           time
> #MaternalAge   14   21
>   #        3 1104  864
>    #       4  875  667
>     #     5   67   53 #less number of observations
>
>
>  BP.stack6 <- BP.stack5[order(BP.stack5$CODEA, BP.stack5$time),]
>  head(BP.stack6)  # very few IDs with  MaternalAge==5
> #       X CODEA Sex MaternalAge Education Birthplace AggScore IntScore
> #1493 3.1     3   2           3         3          1        0        0
> #3202 3.2     3   2           3         3          1        0        0
> #1306 7.1     7   2           4         6          1        0        0
> #3064 7.2     7   2           4         6          1        0        0
> #1    8.1     8   2           4         4          1        0        0
> #2047 8.2     8   2           4         4          1        0        0
>  #         Categ time Obese Overweight hibp
> #1493 Overweight   14     0          0    0
> #3202 Overweight   21     0          1    0
> #1306      Obese   14     0          0    0
> #3064      Obese   21     1          1    0
> #1        Normal   14     0          0    0
> #2047     Normal   21     0          0    0
> BP.stack7<-BP.stack6[BP.stack6$MaternalAge!=5,]
>  BP.stack7$MaternalAge<-factor(as.numeric(as.character(BP.stack7$MaternalAge)
>  
>fit5<-geese(hibp~MaternalAge*time,id=CODEA,data=BP.stack7,family=binomial,corstr="unstructured",scale.fix=TRUE)
> #Error in `contrasts<-`(`*tmp*`, value = contr.funs[1 + isOF[nn]]) :
>  # contrasts can be applied only to factors with 2 or more levels
>
>  with(BP.stack7,table(MaternalAge,time))  #It looks like the combinations
> are still there
> #           time
> #MaternalAge   14   21
>  #         3 1104  864
>    #       4  875  667
>
> It works also with corstr="ar1".   Why do you gave the option
> "unstructured"?
> A.K.
>
>
>
>
>
>
> - Original Message -
> From: rex2013 
> To: r-help@r-project.org
> Cc:
> Sent: Monday, January 7, 2013 6:15 AM
> Subject: Re: [R] random effects model
>
> Hi A.K
>
> Below is the comment I get, not sure why.
>
> BP.sub3 is the stacked data without the missing values.
>
> BP.geese3 <- geese(HiBP~time*MaternalAge,data=BP.sub3,id=CODEA,
> family=binomial, corstr="unstructured", na.action=na.omit)Error in
> `contrasts<-`(`*tmp*`, value = contr.funs[1 + isOF[nn]]) :
>   contrasts can be applied only to factors with 2 or more levels
>
> Even though age has 3 levels; time has 14 years & 21 year

[R] multiple versions of function

2013-01-07 Thread ivo welch
dear R experts:

I want to define a function the calculates the black-scholes value.
it takes 5 named parameters, BS <- function(S,K,dt,rf,sigma) {} .
let's presume I want to be able to call this not only with my 5
numeric vectors BS( sigma=0.3, S=100, K=100, dt=1, rf=0.1 ) and BS(
100, 100, 1, 0.1, 0.3), but also with a data frame that contains the
variables alll in a neat data frame already, BS( data.frame( S=100,
K=100, dt=1, rf=0.1, sigma=0.3 )).  I could of course define BS6 and
BS1, but it would be nice to wrap this functionality into one function
that can do both.

I know that BS has to parse an '...' argument, but there could be a
couple of magical R functions that might make this easier than I would
do it with my planned clunky version. what's the elegant version?

/iaw



Ivo Welch (ivo.we...@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] multiple versions of function

2013-01-07 Thread David Winsemius

On Jan 7, 2013, at 3:57 PM, ivo welch wrote:

> dear R experts:
> 
> I want to define a function the calculates the black-scholes value.
> it takes 5 named parameters, BS <- function(S,K,dt,rf,sigma) {} .
> let's presume I want to be able to call this not only with my 5
> numeric vectors BS( sigma=0.3, S=100, K=100, dt=1, rf=0.1 ) and BS(
> 100, 100, 1, 0.1, 0.3), but also with a data frame that contains the
> variables alll in a neat data frame already, BS( data.frame( S=100,
> K=100, dt=1, rf=0.1, sigma=0.3 )).  I could of course define BS6 and
> BS1, but it would be nice to wrap this functionality into one function
> that can do both.
> 
> I know that BS has to parse an '...' argument, but there could be a
> couple of magical R functions that might make this easier than I would
> do it with my planned clunky version. what's the elegant version?
> 

apply( dfrm, 1, BS)

-- 

David Winsemius
Alameda, 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] plot xaxp issue

2013-01-07 Thread Elaine Kuo
Hello,

Thanks again.
But something wrong with the subset after lm

>
> It is not about drawing plot many times but coding points or graphic objects 
> by some factor. In your case sex.
>
> Instead of
> boyline<-lm(body_weight ~ body_length, boy)
>
> use collective data frame together and subset only one sex.
>
boy<-read.csv("H:/boy_data.csv",header=T)
girl<-read.csv("H:/girl_data.csv",header=T)
together <- rbind(boy, girl)
together$sex <- factor(rep(c("boy", "girl"), c(8,7)))
plot(together$body_length, together$body_weight, ,
col=c("firebrick3","saddlebrown")[as.numeric(together$sex)], )

boyline<-lm(body_weight ~ body_length, data=together, subset(sex=="boy"))
> abline(boyline,col="firebrick3",lwd=2)
>

I updated it aas
>boyline<-lm(body_weight ~ body_length, data=together, 
>subset(together,sex=="boy"))

But an error showed, saying:
error in xj[i] : useless type of  'list'

Please kindly help and thanks again.
BTW, the code below might need update as well.
Elaine


> In case of more than two sexes (some SF authors mentioned it is possible) you 
> can use simple loop.
>
> for (i in 1:2)) {
> subs <- together$sex==levels(together$sex)[i]
> line<-lm(body_weight ~ body_length, data=together, subset=subs)
> abline(line, col=c("firebrick3","saddlebrown")[i],lwd=2)
> }


The original code with advice by Petr
#~~
boy<-read.csv("H:/boy_data.csv",header=T)
girl<-read.csv("H:/girl_data.csv",header=T)
together <- rbind(boy, girl)
together$sex <- factor(rep(c("boy", "girl"), c(8,7)))
plot(together$body_length, together$body_weight, ,
col=c("firebrick3","saddlebrown")[as.numeric(together$sex)], )

boyline<-lm(body_weight ~ body_length, data=together, subset(sex=="boy"))
abline(boyline,col="firebrick3",lwd=2)

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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 xaxp issue

2013-01-07 Thread Elaine Kuo
Hello,

I figured out that the code should be

boyline<-lm(body_weight ~ body_length, data=subset (together,,sex=="boy"))

However, the "" could be omitted if the field name happened to be
numeric, such as 1, 2, or 3.
Please kindly explain why the "" could be omitted for numbers.
Thanks again.

Elaine
On Tue, Jan 8, 2013 at 8:40 AM, Elaine Kuo  wrote:
> boyline<-lm(body_weight ~ body_length, data=together, 
> subset(together,sex=="boy"))

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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 in a abline loop

2013-01-07 Thread Elaine Kuo
Hello

I have data of body length and body weight of people of different skin colors.

I tried to write a code to plot body length and body weight according
to the skin colors.
(Thanks for Petr's advice so far.)

A loop is used but an error shows up in the following code.
It says:
 unexpected '}' in
"
"red3","red3","saddlebrown","coral4","chocolate4","darkblue","navy","grey38")[i],lwd=2)
}"

Please kindly advise how to modify the code.
Thank you.

The code
  data <-read.csv("H:/skincolor.csv",header=T)

  # graph
par(mai=c(1.03,1.03,0.4,0.4))

plot(data$body_weight, data$body_length,
xaxp=c(0,200,4),
yaxp=c(0,200,4),
type="p",
pch=1,lwd=1.0,
cex.lab=1.4, cex.axis=1.2,
font.axis=2,
cex=1.5,
las=1,
bty="l",col=c("yellow","chocolate1","darkorange2",
"red3","saddlebrown","coral4","grey38")[as.numeric(data$skin_color)])



#~~~
##
for (i in 1:7) {
subs <- data$skin_color==levels(data$skin_color)[i]
line<-lm(body_weight~body_length, data=subset(data, subset=subs),
abline(line,col=c("yellow","chocolate1","darkorange2",
"red3","saddlebrown","coral4","grey38")[i],lwd=2)
}

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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 in a abline loop

2013-01-07 Thread Pascal Oettli

Hello,

A ")" is missing.

I guess your code should be

for (i in 1:7){
subs <- data$skin_color==levels(data$skin_color)[i]
line <-lm(body_weight~body_length, data=subset(data, subset=subs)) 
abline(line,col=c("yellow","chocolate1","darkorange2","red3","saddlebrown","coral4","grey38")[i],lwd=2)

}

HTH,
Pascal

Le 08/01/2013 10:23, Elaine Kuo a écrit :

Hello

I have data of body length and body weight of people of different skin colors.

I tried to write a code to plot body length and body weight according
to the skin colors.
(Thanks for Petr's advice so far.)

A loop is used but an error shows up in the following code.
It says:
  unexpected '}' in
"
"red3","red3","saddlebrown","coral4","chocolate4","darkblue","navy","grey38")[i],lwd=2)
 }"

Please kindly advise how to modify the code.
Thank you.

The code
   data <-read.csv("H:/skincolor.csv",header=T)

   # graph
 par(mai=c(1.03,1.03,0.4,0.4))

 plot(data$body_weight, data$body_length,
 xaxp=c(0,200,4),
 yaxp=c(0,200,4),
 type="p",
 pch=1,lwd=1.0,
 cex.lab=1.4, cex.axis=1.2,
 font.axis=2,
 cex=1.5,
 las=1,
 bty="l",col=c("yellow","chocolate1","darkorange2",
"red3","saddlebrown","coral4","grey38")[as.numeric(data$skin_color)])


 
#~~~
 ##
 for (i in 1:7) {
 subs <- data$skin_color==levels(data$skin_color)[i]
 line<-lm(body_weight~body_length, data=subset(data, subset=subs),
 abline(line,col=c("yellow","chocolate1","darkorange2",
"red3","saddlebrown","coral4","grey38")[i],lwd=2)
 }

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] error in a abline loop

2013-01-07 Thread Elaine Kuo
Hello

Thanks.
It is right a ) is missing for lm.
With the modification, I ran the code but the command seems not
closed, with + instead of >.
Please kindly help.

The code modified
for (i in 1:7) {
subs <- data$skin_color==levels(data$skin_color)[i]
line<-lm(body_weight~body_length, data=subset(data, subset=subs))
abline(line,col=c("yellow","chocolate1","darkorange2",
"red3","saddlebrown","coral4","grey38")[i],lwd=2) )
}

Elaine



On Tue, Jan 8, 2013 at 10:03 AM, arun  wrote:
>
> HI,
>
> A possible guess ( with no data):
> for (i in 1:7) {
> subs <- data$skin_color==levels(data$skin_color)[i]
> line<-lm(body_weight~body_length, data=subset(data, subset=subs)) 
> #closing parentheses
> abline(line,col=c("yellow","chocolate1","darkorange2",
> "red3","saddlebrown","coral4","grey38")[i],lwd=2)
> }
>
>
> A.K.
>
>
>
> - Original Message -
> From: Elaine Kuo 
> To: r-help@r-project.org
> Cc:
> Sent: Monday, January 7, 2013 8:23 PM
> Subject: [R] error in a abline loop
>
> Hello
>
> I have data of body length and body weight of people of different skin colors.
>
> I tried to write a code to plot body length and body weight according
> to the skin colors.
> (Thanks for Petr's advice so far.)
>
> A loop is used but an error shows up in the following code.
> It says:
> unexpected '}' in
> "
> "red3","red3","saddlebrown","coral4","chocolate4","darkblue","navy","grey38")[i],lwd=2)
> }"
>
> Please kindly advise how to modify the code.
> Thank you.
>
> The code
>   data <-read.csv("H:/skincolor.csv",header=T)
>
>   # graph
> par(mai=c(1.03,1.03,0.4,0.4))
>
> plot(data$body_weight, data$body_length,
> xaxp=c(0,200,4),
> yaxp=c(0,200,4),
> type="p",
> pch=1,lwd=1.0,
> cex.lab=1.4, cex.axis=1.2,
> font.axis=2,
> cex=1.5,
> las=1,
> bty="l",col=c("yellow","chocolate1","darkorange2",
> "red3","saddlebrown","coral4","grey38")[as.numeric(data$skin_color)])
>
>
> 
> #~~~
> ##
> for (i in 1:7) {
> subs <- data$skin_color==levels(data$skin_color)[i]
> line<-lm(body_weight~body_length, data=subset(data, subset=subs),
> abline(line,col=c("yellow","chocolate1","darkorange2",
> "red3","saddlebrown","coral4","grey38")[i],lwd=2)
> }
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/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] error in a abline loop

2013-01-07 Thread Pascal Oettli

Hi,

Now, you have too many ")" in the code. The last line should be (as 
written in the code I sent you):


> abline(line,col=c("yellow","chocolate1","darkorange2",
> "red3","saddlebrown","coral4","grey38")[i],lwd=2)

HTH,
Pascal


Le 08/01/2013 11:24, Elaine Kuo a écrit :

Hello

Thanks.
It is right a ) is missing for lm.
With the modification, I ran the code but the command seems not
closed, with + instead of >.
Please kindly help.

The code modified
for (i in 1:7) {
 subs <- data$skin_color==levels(data$skin_color)[i]
 line<-lm(body_weight~body_length, data=subset(data, subset=subs))
 abline(line,col=c("yellow","chocolate1","darkorange2",
"red3","saddlebrown","coral4","grey38")[i],lwd=2) )
 }

Elaine



On Tue, Jan 8, 2013 at 10:03 AM, arun  wrote:


HI,

A possible guess ( with no data):
for (i in 1:7) {
 subs <- data$skin_color==levels(data$skin_color)[i]
 line<-lm(body_weight~body_length, data=subset(data, subset=subs)) #closing 
parentheses
 abline(line,col=c("yellow","chocolate1","darkorange2",
"red3","saddlebrown","coral4","grey38")[i],lwd=2)
 }


A.K.



- Original Message -
From: Elaine Kuo 
To: r-help@r-project.org
Cc:
Sent: Monday, January 7, 2013 8:23 PM
Subject: [R] error in a abline loop

Hello

I have data of body length and body weight of people of different skin colors.

I tried to write a code to plot body length and body weight according
to the skin colors.
(Thanks for Petr's advice so far.)

A loop is used but an error shows up in the following code.
It says:
unexpected '}' in
"
"red3","red3","saddlebrown","coral4","chocolate4","darkblue","navy","grey38")[i],lwd=2)
 }"

Please kindly advise how to modify the code.
Thank you.

The code
   data <-read.csv("H:/skincolor.csv",header=T)

   # graph
 par(mai=c(1.03,1.03,0.4,0.4))

 plot(data$body_weight, data$body_length,
 xaxp=c(0,200,4),
 yaxp=c(0,200,4),
 type="p",
 pch=1,lwd=1.0,
 cex.lab=1.4, cex.axis=1.2,
 font.axis=2,
 cex=1.5,
 las=1,
 bty="l",col=c("yellow","chocolate1","darkorange2",
"red3","saddlebrown","coral4","grey38")[as.numeric(data$skin_color)])


 
#~~~
 ##
 for (i in 1:7) {
 subs <- data$skin_color==levels(data$skin_color)[i]
 line<-lm(body_weight~body_length, data=subset(data, subset=subs),
 abline(line,col=c("yellow","chocolate1","darkorange2",
"red3","saddlebrown","coral4","grey38")[i],lwd=2)
 }

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] multiple versions of function

2013-01-07 Thread ivo welch
hi david---can you give just a little more of an example?  the
function should work with call by order, call by name, and data frame
whose columns are the names.  /iaw



Ivo Welch (ivo.we...@gmail.com)


On Mon, Jan 7, 2013 at 4:25 PM, David Winsemius  wrote:
>
> On Jan 7, 2013, at 3:57 PM, ivo welch wrote:
>
>> dear R experts:
>>
>> I want to define a function the calculates the black-scholes value.
>> it takes 5 named parameters, BS <- function(S,K,dt,rf,sigma) {} .
>> let's presume I want to be able to call this not only with my 5
>> numeric vectors BS( sigma=0.3, S=100, K=100, dt=1, rf=0.1 ) and BS(
>> 100, 100, 1, 0.1, 0.3), but also with a data frame that contains the
>> variables alll in a neat data frame already, BS( data.frame( S=100,
>> K=100, dt=1, rf=0.1, sigma=0.3 )).  I could of course define BS6 and
>> BS1, but it would be nice to wrap this functionality into one function
>> that can do both.
>>
>> I know that BS has to parse an '...' argument, but there could be a
>> couple of magical R functions that might make this easier than I would
>> do it with my planned clunky version. what's the elegant version?
>>
>
> apply( dfrm, 1, BS)
>
> --
>
> David Winsemius
> Alameda, 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] error in a abline loop

2013-01-07 Thread arun

HI,

A possible guess ( with no data):
for (i in 1:7) {
    subs <- data$skin_color==levels(data$skin_color)[i]
    line<-lm(body_weight~body_length, data=subset(data, subset=subs)) #closing 
parentheses
    abline(line,col=c("yellow","chocolate1","darkorange2",
"red3","saddlebrown","coral4","grey38")[i],lwd=2)  
    }


A.K.



- Original Message -
From: Elaine Kuo 
To: r-help@r-project.org
Cc: 
Sent: Monday, January 7, 2013 8:23 PM
Subject: [R] error in a abline loop

Hello

I have data of body length and body weight of people of different skin colors.

I tried to write a code to plot body length and body weight according
to the skin colors.
(Thanks for Petr's advice so far.)

A loop is used but an error shows up in the following code.
It says:
unexpected '}' in
"    
"red3","red3","saddlebrown","coral4","chocolate4","darkblue","navy","grey38")[i],lwd=2)
    }"

Please kindly advise how to modify the code.
Thank you.

The code
  data <-read.csv("H:/skincolor.csv",header=T)

  # graph
    par(mai=c(1.03,1.03,0.4,0.4))

    plot(data$body_weight, data$body_length,
    xaxp=c(0,200,4),
    yaxp=c(0,200,4),
    type="p",
    pch=1,lwd=1.0,
    cex.lab=1.4, cex.axis=1.2,
    font.axis=2,
    cex=1.5,
    las=1,
    bty="l",col=c("yellow","chocolate1","darkorange2",
"red3","saddlebrown","coral4","grey38")[as.numeric(data$skin_color)])


    
#~~~
    ##
    for (i in 1:7) {
    subs <- data$skin_color==levels(data$skin_color)[i]
    line<-lm(body_weight~body_length, data=subset(data, subset=subs),
    abline(line,col=c("yellow","chocolate1","darkorange2",
"red3","saddlebrown","coral4","grey38")[i],lwd=2)
    }

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Conditional Statistics

2013-01-07 Thread Joseph Norman Thomson
Hello,

I am a new user of R. I am coming from SAS and do statistics on stock
market data, economic data, and social data. My question is this: How
can you get the mean, standard dev, etc. of a variable based on a
conditional statement on either the same variable or a different
variable in the same data set? So if I had the closing prices of the
S&P from 01/01/1990-12/31/1990, how could I get the average price of
the S&P from 02/01/1990-03/15/1990? Or the average price of the S&P on
Mondays (assuming a dummy var is created for 1 = Monday, 0 = else). I
understand that you can create subsets and new data sets based on the
conditional statements; but is there an easier way to do this by
typing a line into the mean() statement? That was extremely easy in
SAS where you could say:

proc means data=sp500;
var price;
where monday = 1;

Thank you for your help.

Joe

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Conditional Statistics

2013-01-07 Thread Simon Blomberg
You can use the tapply function to do this. You can't type a line into 
the mean statement. (See ?mean for what you can type in there). The 
general approach is to have a vector of data (stock prices) and a 
categorical variable (day of week). Then break up the data vector 
according to the levels in the categorical variable, and calculate the 
mean values:


Weekmeans <- tapply(data.vector, catvariable, mean)

This will give you the means for all days. If you really just want one 
mean (just monday), you could do:


Monmean <- mean(data.vector[catvariable=="Monday"])

Similarly, if you want the standard deviation for each day of the week, 
you would use:


WeekSD <- tapply(data.vector, catvariable, sd)
MonSD <- sd(data.vector[catvariable=="Monday"])

You will find that some things that are easy in SAS require a little 
more thought in R, and vice versa. Certainly, the philosophical approach 
to data analysis in R is different to that in SAS. There are a couple of 
books for R for SAS users. They might help you.


Cheers,

Simon.

On 08/01/13 11:17, Joseph Norman Thomson wrote:

Hello,

I am a new user of R. I am coming from SAS and do statistics on stock
market data, economic data, and social data. My question is this: How
can you get the mean, standard dev, etc. of a variable based on a
conditional statement on either the same variable or a different
variable in the same data set? So if I had the closing prices of the
S&P from 01/01/1990-12/31/1990, how could I get the average price of
the S&P from 02/01/1990-03/15/1990? Or the average price of the S&P on
Mondays (assuming a dummy var is created for 1 = Monday, 0 = else). I
understand that you can create subsets and new data sets based on the
conditional statements; but is there an easier way to do this by
typing a line into the mean() statement? That was extremely easy in
SAS where you could say:

proc means data=sp500;
var price;
where monday = 1;

Thank you for your help.

Joe

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


--
Simon Blomberg, BSc (Hons), PhD, MAppStat, AStat.
Lecturer and Consultant Statistician
School of Biological Sciences
The University of Queensland
St. Lucia Queensland 4072
Australia
T: +61 7 3365 2506
email: S.Blomberg1_at_uq.edu.au
http://www.evolutionarystatistics.org

Policies:
1.  I will NOT analyse your data for you.
2.  Your deadline is your problem.

Statistics is the grammar of science - Karl Pearson.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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 xaxp issue

2013-01-07 Thread PIKAL Petr
Hi

I actually do not know what you are talking about. I believe that the code I 
provided works (if i did not make typing mistake), but you failed to provide 
data to test so I made my own data and subset inside lm works as expected.

fake<-data.frame(weight<-rnorm(15), length=weight*3+rnorm(15), 
sex=rep(c("b","g"), c(7,8)))
names(fake)[1]<-"weight"
plot(fake$length, fake$weight, col=as.numeric(fake$sex))
bl<-lm(weight~length, data=fake, subset=sex=="b")
abline(bl)
gl<-lm(weight~length, data=fake, subset=sex=="g")
abline(gl, col=2)

Beware that R is smart but not smart enough to correct user typing mistakes 
(missing parntheses, missing commas, etc.). Typing mistakes are probably the 
most common source of failing code.

Regards

Petr

> -Original Message-
> From: Elaine Kuo [mailto:elaine.kuo...@gmail.com]
> Sent: Tuesday, January 08, 2013 2:03 AM
> To: PIKAL Petr
> Cc: r-help@r-project.org
> Subject: Re: [R] plot xaxp issue
> 
> Hello,
> 
> I figured out that the code should be
> 
> boyline<-lm(body_weight ~ body_length, data=subset
> (together,,sex=="boy"))
> 
> However, the "" could be omitted if the field name happened to be
> numeric, such as 1, 2, or 3.
> Please kindly explain why the "" could be omitted for numbers.
> Thanks again.
> 
> Elaine
> On Tue, Jan 8, 2013 at 8:40 AM, Elaine Kuo 
> wrote:
> > boyline<-lm(body_weight ~ body_length, data=together,
> > subset(together,sex=="boy"))

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Have problem to do loop to generate transformed chi-squared variates

2013-01-07 Thread Rolf Turner



Homework?  We don't do people's homework for them.

cheers,

Rolf Turner

On 01/08/2013 05:13 AM, Agnes Ayang wrote:

Hello R-helpers,

I need to generate standard variates normal to 'create' chi-squared variates. 
To make you more understand,

(1)   a<-rnorm(3,0,1)

*after do (1), I need to squared and summed the three values. My problem is, 
how am I going to continue the programming if I had to repeat the process for 
15 times, which in the end I will get 15 values from the whole programme.Hope 
you can help.


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