[R] Repeated Measures design using lme
Hi, I have what I believe is a repeated-measures dataset that I'm trying to analyze using lme(). This is *not* homework, but an exercise in my trying to self-teach myself repeated-measure ANOVA for other *real* datasets that I have and that are extremely similar to the following design. I'm fairly sure the dataset described below would work with lme() -- but it'd be great if anybody can confirm that after I describe the dataset below) The study involves measuring the effect of a drug on blood pressure. There were 16 patients in all and 6 replicate measures per patient of their blood pressure on one week (one measure per day). Two weeks later, a drug was introduced to 8 randomly selected patients in such a way that I had equal representation of the 4 age groups among the two treatment groups. Then, another two weeks later, 6 replicate measures per patient (per day) of blood pressure was retaken. So each patient had 12 total measures whether they were in the treatment group or in the control group (6 reads (R1-R6) in the baseline-week and another 6 reads (R1-R6) in the post-treatment week). So, Background: 16 patients Response measure: Blood pressure Fixed Factor: 4 Age groups Fixed Factor: Drug vs. NoDrug Random factor: Day of the read (i.e. 6 replicate reads (R1-R6) at the baseline time, and 6 replicate reads (R1-R6) after the drug has had time to take effect) Random Factor: Subjects 1-16 PatientAgeGroupBP(Blood Pressure)Read (replicate reads) Pre/PostTreatmentWeekGroup 120-2983 R1 preTreat 220-2981 R1 preControl 320-2974 R1 preTreat 420-2985 R1 preControl 530-3982 R1 preTreat ……………… 320-2974 R2 preTreat ……………… 120-2983 R1 post Treat 220-2982 R1 post Control 320-2986 R1 post Treat 420-2984 R1 post Control ……………… I'm trying to do an analysis of variance to decide whether there is a measurable change in blood pressure between the Treat and Control groups. Another issue is that some of the 16 patients didn't get all 6 replicate reads in their pre/post treatment weeks, so I need to include the na.omit function. What I think I'm having the most trouble with is the repeated reads (R1 through R6) in the pre/post treatment weeks. I'm fairly sure this is a random variable -- their order or identify (R1 in pre-treatment week has no relation to R1 in the post-treatment week, etc). By placing Read as a random variable, am I covering myself there? If I execute: summary(lme(BP ~ Group, random = ~ 1 | Patient, data = bloodpress, na.action=na.omit)) I get a result, but I'm not sure it's correct -- do I need to tell the model about the Read factor (R1-R6 in pre/post weeks)? I'm really trying to set the right form of the lme() function call to decide 1) if there is a statistical difference between the Treat/Control groups and, 2) if one takes into account AgeGroup, is there a statistical difference between Treat/Control Groups, and finally 3) if I don't see a statistical difference, can someone recommend an R function that might solve the supplemental question, given the noise in day-to-day blood pressure reads, and given that I wanted to have enough statistical power to observe a say, a 5% benefit in blood pressure, how many additional reads or patients I would need. Basically, is lme() the proper function, and can someone offer any pointers on what my call to this function should look like to make the above to determinations? Thanks! -Scott __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE
[R] Parsing results from boot
This probably has a super easy answer...but I claim newbie status! (I did search help lists but this question is hard to isolate keyword-wise) Basically, I'm trying to figure out how to parse the results from executing boot(). I'm mainly interested in assigning the standard error estimate to a scalar variable. For example: -- b-c(100,100,120,130,1000,1200,1100,1150,125) b.boot - boot(b, function(b,i) median(b[i]), R=1000) b.boot ORDINARY NONPARAMETRIC BOOTSTRAP Call: boot(data = b, statistic = function(b, i) median(b[i]), R = 1000) Bootstrap Statistics : original biasstd. error t1* 130 354.75450.9763 -- I'm interested in the value for std.error (i.e. 450.9763). Now executing the folling: -- summary(b.boot) Length Class Mode t0 1 -none- numeric t 1000 -none- numeric R1 -none- numeric data 9 -none- numeric seed 626 -none- numeric statistic1 -none- function sim 1 -none- character call 4 -none- call stype1 -none- character strata 9 -none- numeric weights 9 -none- numeric --- seems to imply that it's not actually assigned to an output variable...but it is shown when I just type b.boot Where is the std.error being assigned? and how to I access it? Thanks in advance... -Scott __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Ordering bars in barplots
I have a dataframe such that when I enter the dataframe name at the R prompt and see the data, the order of the data is correct (ie. what I want - ordered numerically by the factor, class). The table is akin to the following: df: (dataframe = df) uniqueID class age a 132 b 121 c 113 d 111 e 315 f 316 g 331 h 325 i 423 j 432 k 431 l 711 m 76 n 720 o 1019 p 1025 q 1042 . . . uniqueID and class are factors. Yesterday I think I managed to figure out how to arrange the column class in numerical order (when it lists) even though it's a character factor. (although feel free to pipe in an answer on that too - since I was trying a lot of different things!) Now, when I execute the following command: barplot(by(df[,age],df$class,FUN=function(x) quantile(x,0.95))) the x-labels or class of each bar is in *character* order, not *numerical order. Is there a way to rearrage this so that each bar label (categorical x-label) is in numerical order? Thanks in advance!!! Best, Scott [[alternative HTML version deleted]] __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Collapsing a factor in R
I'm trying to collapse the following table along the sub-group factor. In this case, collapsing means taking the average ages within a subgroup and creating a new table. I seem to be running into trouble trying to create this new data frame. I can use the ave() function to find averages within a subgroup but how do I maintain the Group/Gender factors after collapsing? (see bottom table) Can anyone help?? Thanks in advance!! Here's the starting point (table) Age Group SubGroup Gender 12 g1 a f 32 g1 a f 81 g1 a f 63 g1 b m 24 g1 b m 24 g1 b m 70 g1 c f 82 g1 c f 71 g1 c f 61 g1 c f 25 g2 d m 29 g2 d m 43 g2 e f 17 g2 e f 42 g2 e f 55 g2 e f 52 g2 f f 19 g2 f f 12 g2 f f The following is what I'm trying to achieve (as a new dataframe) AvgAgeSubgroup Group SubgroupGender 41.67 g1 a f 37 g1 b m 71 g1 c f 27 g2 d m 39.25 g2 e f 27.67 g2 f f Can anyone help? -Scott [[alternative HTML version deleted]] __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] How to create unique factor from two factors? + Boostrap Q
This might be easy but I'm very new to R and this question doesn't seem to have any nice keywords that bring up relevant search results when I search the CRAN search engine. Therefore, I'll plead (as I have in the recent past) Newbie status. I have a data frame with two factors (Factor 1 and 2) which together specify another unique level. I want to create a third factor in the data frame that captures this uniqueness. For example, say I had dataframe, Df, with Factors, 1 and 2. I want to create Factor 3 and add it to my Df dataframe. i.e. Df dataframe: WANT TO Row# Factor1 Factor2 CREATE THIS: Factor 3Data 11 11 23 21 22 43 31 22 19 41 22 11 51 43 3 61 43 13 73 14 52 83 14 12 93 14 9 10 3 35 21 11 3 35 43 12 4 16 32 13 4 16 18 14 4 27 52 15 4 27 21 and of course, I'm trying to create Factor 3 without loops.. My end goal here (which I add because maybe I don't need to create Factor 3 (although I'm still curious)), is to bootstrap sample Factor 3. I want to repeatedly grab, say, 3 levels of Factor 3, and take the mean of those levels (e.g. say in my first bootstrap sample, I grab levels 2,4, and 7 from Factor 3, then I want to take the mean of rows, 2,3,4,7,8,9,14,15). Of course, each sample from Factor 3 for my bootstrap will most likely have a differing number of rows since my experiment is not balanced. I'm not sure if this is an issue yet when I try to implement the boot function in R (I haven't gotten to that point yet). Any help is appreciated! Thanks! -Scott Scott Norton, Ph.D. Engineering Manager Nanoplex Technologies, Inc. 2375 Garcia Ave. Mountain View, CA 94043 www.nanoplextech.com http://www.nanoplextech.com/ [EMAIL PROTECTED] [[alternative HTML version deleted]] __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
RE: [R] outer function problems
Thanks Spencer and Tom for your help!
Besides the other errors, I realized last night that I'm making a
fundmental error in my interpretation of the outer function. The following
short code snippet highlights my confusion.
f-function(A,B) { sum(A+B) }
outer(1:3,2:4,f)
[,1] [,2] [,3]
[1,] 45 45 45
[2,] 45 45 45
[3,] 45 45 45
I had *thought* that outer() would give:
[,1] [,2] [,3]
[1,] 3 45
[2,] 4 56
[3,] 5 67
ie. take each combination from A = 1,2,3; B=2,3,4 such as A=1,B=2 put it in
the sum function, get [1,1]=3 ...
Then grab A[2]=2,B[1]=2, put them in the sum() function to get [2,1]=4,
etc... That seems to be the way the instructions explain outer, i.e.
element-by-element computation of FUN()
Description:
The outer product of the arrays 'X' and 'Y' is the array 'A' with
dimension 'c(dim(X), dim(Y))' where element 'A[c(arrayindex.x,
arrayindex.y)] = FUN(X[arrayindex.x], Y[arrayindex.y], ...)'.
Since my interpretation is *definitely* wrong, could someone put in words
how OUTER handles the argument vectors and the functional call with
reference to the preceding example?
Also, what I need to happen in my code is to actually take each combination
of elements from vectors, A and B, and feed them repeatedly into a
function, generating a matrix of results. How then do I do that?
Thanks in advance!!!
-Scott
Scott Norton, Ph.D.
Engineering Manager
Nanoplex Technologies, Inc.
2375 Garcia Ave.
Mountain View, CA 94043
www.nanoplextech.com
-Original Message-
From: Spencer Graves [mailto:[EMAIL PROTECTED]
Sent: Tuesday, October 28, 2003 8:13 AM
To: Scott Norton
Cc: [EMAIL PROTECTED]
Subject: Re: [R] outer function problems
I don't know that this is your problem, but I see a potential
scoping issue: It is not obvious to me where Dk is getting n0 and w.
I've solved this kind of problem in the past by declaring n0 and w as
explicit arguments to Dk and then passing them explicitly via ... in
outer. In general, I prefer to avoid accessing globals from within
functions. This may not help you here, but it might help in the future.
hope this helps. spencer graves
Scott Norton wrote:
I'm pulling my hair (and there's not much left!) on this one. Basically I'm
not getting the same result t when I step through the program and
evaluate
each element separately than when I use the outer() function in the
FindLikelihood() function below.
Here's the functions:
Dk- function(xk,A,B)
{
n0 *(A*exp(-0.5*(xk/w)^2) + B)
}
FindLikelihood - function(Nk)
{
A - seq(0.2,3,by=0.2)
B - seq(0.2,3,by=0.2)
k -7
L - outer(A, B, function(A,B) sum( (Nk*log(Dk(seq(-k,k),A,B))) -
Dk(seq(-k,k),A,B) ))
return(L)
}
where Nk - c(70 , 67 , 75 , 77 , 74 ,102, 75, 104 , 94 , 74 , 78 , 79 ,
83
, 73 , 76)
Here's an excerpt from my debug session..
Nk
[1] 70 67 75 77 74 102 75 104 94 74 78 79 83 73 76
debug(FindLikelihood)
L-FindLikelihood(Nk)
debugging in: FindLikelihood(Nk)
debug: {
A - seq(0.2, 3, by = 0.2)
B - seq(0.2, 3, by = 0.2)
k - 7
L - outer(A, B, function(A, B) sum((Nk * log(Dk(seq(-k,
k), A, B))) - Dk(seq(-k, k), A, B)))
return(L)
}
Browse[1] n
debug: A - seq(0.2, 3, by = 0.2)
Browse[1] n
debug: B - seq(0.2, 3, by = 0.2)
Browse[1] n
debug: k - 7
Browse[1] n
debug: L - outer(A, B, function(A, B) sum((Nk * log(Dk(seq(-k, k),
A, B))) - Dk(seq(-k, k), A, B)))
Browse[1] sum((Nk * log(Dk(seq(-k, k),0.2,0.2))) - Dk(seq(-k, k), 0.2,
0.2)) # WHY DOES THIS LINE GIVE ME THE CORRECT RESULT WHEN I
SUBSTITUTE
0.2, 0.2 FOR A AND B
[1] 2495.242
Browse[1] outer(A, B, function(A, B) sum((Nk * log(Dk(seq(-k, k),
+ A, B))) - Dk(seq(-k, k), A, B)))
[,1] [,2] [,3] [,4] [,5] [,6] [,7]
[,8]
[1,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48# BUT ELEMENT (1,1) WHICH SHOULD ALSO BE (A,B) = (0.2, 0.2),
GIVES THE INCORRECT RESULT
[2,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[3,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[4,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[5,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[6,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[7,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[8,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[9,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[10,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[11,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[12,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[13,] 58389.48 58389.48 58389.48 58389.48 58389.48
[R] outer function problems
I'm pulling my hair (and there's not much left!) on this one. Basically I'm
not getting the same result t when I step through the program and evaluate
each element separately than when I use the outer() function in the
FindLikelihood() function below.
Here's the functions:
Dk- function(xk,A,B)
{
n0 *(A*exp(-0.5*(xk/w)^2) + B)
}
FindLikelihood - function(Nk)
{
A - seq(0.2,3,by=0.2)
B - seq(0.2,3,by=0.2)
k -7
L - outer(A, B, function(A,B) sum( (Nk*log(Dk(seq(-k,k),A,B))) -
Dk(seq(-k,k),A,B) ))
return(L)
}
where Nk - c(70 , 67 , 75 , 77 , 74 ,102, 75, 104 , 94 , 74 , 78 , 79 , 83
, 73 , 76)
Here's an excerpt from my debug session..
Nk
[1] 70 67 75 77 74 102 75 104 94 74 78 79 83 73 76
debug(FindLikelihood)
L-FindLikelihood(Nk)
debugging in: FindLikelihood(Nk)
debug: {
A - seq(0.2, 3, by = 0.2)
B - seq(0.2, 3, by = 0.2)
k - 7
L - outer(A, B, function(A, B) sum((Nk * log(Dk(seq(-k,
k), A, B))) - Dk(seq(-k, k), A, B)))
return(L)
}
Browse[1] n
debug: A - seq(0.2, 3, by = 0.2)
Browse[1] n
debug: B - seq(0.2, 3, by = 0.2)
Browse[1] n
debug: k - 7
Browse[1] n
debug: L - outer(A, B, function(A, B) sum((Nk * log(Dk(seq(-k, k),
A, B))) - Dk(seq(-k, k), A, B)))
Browse[1] sum((Nk * log(Dk(seq(-k, k),0.2,0.2))) - Dk(seq(-k, k), 0.2,
0.2)) # WHY DOES THIS LINE GIVE ME THE CORRECT RESULT WHEN I SUBSTITUTE
0.2, 0.2 FOR A AND B
[1] 2495.242
Browse[1] outer(A, B, function(A, B) sum((Nk * log(Dk(seq(-k, k),
+ A, B))) - Dk(seq(-k, k), A, B)))
[,1] [,2] [,3] [,4] [,5] [,6] [,7]
[,8]
[1,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48# BUT ELEMENT (1,1) WHICH SHOULD ALSO BE (A,B) = (0.2, 0.2),
GIVES THE INCORRECT RESULT
[2,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[3,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[4,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[5,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[6,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[7,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[8,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[9,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[10,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[11,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[12,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[13,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[14,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[15,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
58389.48
[,9][,10][,11][,12][,13][,14][,15]
[1,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
[2,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
[3,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
[4,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
[5,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
[6,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
[7,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
[8,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
[9,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
[10,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
[11,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
[12,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
[13,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
[14,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
[15,] 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48 58389.48
Browse[1]
As commented above, when I evaluate a single A,B element (i.e. A=0.2,
B=0.2) I get a different result than when I use OUTER() which should also be
evaluating at A=0.2, B=0.2??
Any help appreciated. I know I'm probably doing something overlooking
something simple, but can anyone point it out???
Thanks!
-Scott
Scott Norton, Ph.D.
Engineering Manager
Nanoplex Technologies, Inc.
2375 Garcia Ave.
Mountain View, CA 94043
www.nanoplextech.com
[[alternative HTML version deleted]]
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
[R] returning dynamic variable names from function
Within a function I'm assigning dynamic variable names and values to them
using the assign function. I want to pass back the results but am
uncertain how to do this.
Basically, my function reads a number of data files and uses the filename of
each file as the variable name for a list-to-become-dataframe. I want then
to pass all these lists back, but again, the names of the lists to pass back
have been assign-ed to the filename. See code snippet below, especially
the *** WHAT DO I PUT HERE *** part!
Thanks in advance for any help!!!
-Scott
CODE SNIPPET BELOW:
loadFiles - function()
{
fullFileNames - choose.files(filters = c(Fluor file (*.data,*.Dat),
*.dat;*.Dat))
numFiles - length(fullFileNames)
fileNames - basename(fullFileNames) # removes the all the paths from the
full filenames
splitNames - strsplit(fileNames,_)# create a list of strings
separated by the _ character
dfNames - sapply(splitNames,[,1)
# or I could use sapply(splitNames, function (x) x[1]) OR EVEN
unlist(lapply(splitNames, function(x)x[1])) -- [Thanks to Andy Liaw,
Simon Blomberg, Gabor Grothendieck,
#James Holtman and Robert Keefe for their helpful responses on this
question - See posting on returning the first elements of a list of vectors]
#
for (i in 1:numFiles)
{
assign(dfNames[[i]],ReadFileCreateDF(fullFileNames[[i]]))
# we use assign since we want the quoted string in dfNames to be assigned
the value of the function call
# the ReadFileCreateList is a function that reads an
individual file and sends back a number of parameters in a list
}
return(list( *** WHAT DO I PUT HERE***))
}
Scott Norton, Ph.D.
Engineering Manager
Nanoplex Technologies, Inc.
2375 Garcia Ave.
Mountain View, CA 94043
www.nanoplextech.com
[[alternative HTML version deleted]]
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
[R] indexing a particular element in a list of vectors
I have a list of character vectors. I'm trying to see if there is a way (in a single line, without a loop) to pull out the first element of all the vectors contained in the list. listOfVectors[1:length(listOfVectors][1] doesn't work. == If you want more details.. Here is my listOfVectors which is called uuu uuu[order(uuu)] [[1]] [1] pt1pg multi.expr [[2]] [1] 1ngml fluor.expr [[3]] [1] 1pgml fluor.expr [[4]] [1] 10ng ml fluor.expr [[5]] [1] 10pg ml fluor.expr I'm basically interested in getting the following elements: pt1pg, 1ng, 1pg, etc. as a list (or character vector) Note, this might have an obvious solution but I claim Newbie status! Thanks in advance! -Scott Scott Norton, Ph.D. Engineering Manager Nanoplex Technologies, Inc. 2375 Garcia Ave. Mountain View, CA 94043 www.nanoplextech.com [[alternative HTML version deleted]] __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
