Re: [R] x,y,z table to matrix with x as rows and y as columns
jiho jo.irisson at gmail.com writes: Hello all, I am sure I am missing something obvious but I cannot find the function I am looking for. I have a data frame with three columns: X, Y and Z, with X and Y being grid coordinates and Z the value associated with these coordinates. I want to transform this data frame in a matrix of Z values, on the grid defined by X and Y (and, as a plus, fill the X.Y combinations which do no exist in the original data frame with NAs in the resulting matrix). I could do this manually but I guess the appropriate function should be somewhere around. I just can't find it. I have used xtabs in the past, but xtabs returns 0 instead of NA, which makes great for cross-tabulation, the offending line is x[is.na(x)] - 0. So you would need to either modify the xtabs function or trust that z is never 0 and replace all 0's with NA. mydat - expand.grid(x=1:5, y=1:5) mydat - data.frame(mydat, z=rnorm(25)) mydat$z[sample(1:25,4)] - NA mytab - xtabs(z~x+y, mydat) Mark Lyman __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] x,y,z table to matrix with x as rows and y as columns
I am sure I am missing something obvious but I cannot find the
function I am looking for. I have a data frame with three columns: X,
Y and Z, with X and Y being grid coordinates and Z the value
associated with these coordinates. I want to transform this data
frame in a matrix of Z values, on the grid defined by X and Y (and,
as a plus, fill the X.Y combinations which do no exist in the
original data frame with NAs in the resulting matrix). I could do
this manually but I guess the appropriate function should be
somewhere around. I just can't find it.
Immediately after my last post I realized there was a much better solution
mydat - expand.grid(x=1:5, y=1:5)
mydat - data.frame(mydat, z=rnorm(25))
mydat$z[sample(1:25,4)] - NA
data2mat - function(x, y, z)
+ {
+ out - matrix(unlist(split(z, interaction(x,y))), ncol=length(unique(y)))
+ dimnames(out) - list(unique(x), unique(y))
+ out
+ }
with(mydat, data2mat(x, y, z))
Mark Lyman
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] problems with character objects and calls to list()
Really I'd like the call to list() to behave as though the text had been entered directly so that you get list(1:2, 3:4, 5:6) [[1]] [1] 1 2 [[2]] [1] 3 4 [[3]] [1] 5 6 eval(parse(text=paste(list(,to.convert,),sep=))) [[1]] [1] 1 2 [[2]] [1] 3 4 [[3]] [1] 5 6 [[4]] [1] 7 8 [[5]] [1] 9 10 [[6]] [1] 11 12 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Adding points to a wireframe with conditioning variable
I would like to add points to a wireframe but with a conditioning variable. I found a solution for this without a conditioning variable here, http://finzi.psych.upenn.edu/R/Rhelp02a/archive/65321.html. Does anyone know how to plot a wireframe conditioned on a variable and add the points conditioned on the same variable, similar to the solution at the link above? __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Adding points to a wireframe with conditioning variable
Deepayan Sarkar deepayan.sarkar at gmail.com writes:
On 7/5/07, Mark Lyman mark.lyman at gmail.com wrote:
I would like to add points to a wireframe but with a conditioning
variable. I
found a solution for this without a conditioning variable here,
http://finzi.psych.upenn.edu/R/Rhelp02a/archive/65321.html. Does anyone
know
how to plot a wireframe conditioned on a variable and add the points
conditioned on the same variable, similar to the solution at the link
above?
Depends on what form you have your points in. Inside a panel function,
packet.number() will give you the packet number in sequential order
(column major order if you think of the conditioning variables as
defining a multiway cross-tabulation), and which.packet() will give
you a vector with the current level of each conditioning variable. You
can use these to extract the appropriate subset of points.
-Deepayan
Thank you Deepayan. I modified the panel function you wrote in the post
referenced above according to your suggestion. Here is the panel function, I
came up with
panel.3dwire.points - function(x, y, z, xlim, ylim, zlim, xlim.scaled,
ylim.scaled, zlim.scaled, pts, ...)
{
panel.3dwire(x=x, y=y, z=z, xlim=xlim, ylim=ylim, zlim=zlim,
xlim.scaled=xlim.scaled, ylim.scaled=ylim.scaled,
zlim.scaled=zlim.scaled, ...)
pts.sub - subset(pts, g==levels(g)[which.packet()])
xx - xlim.scaled[1] + diff(xlim.scaled)*(pts.sub$x - xlim[1])/diff
(xlim)
yy - ylim.scaled[1] + diff(ylim.scaled)*(pts.sub$y - ylim[1])/diff
(ylim)
zz - zlim.scaled[1] + diff(zlim.scaled)*(pts.sub$z - zlim[1])/diff
(zlim)
panel.3dscatter(x=xx, y=yy, z=zz, xlim=xlim, ylim=ylim, zlim=zlim,
xlim.scaled=xlim.scaled, ylim.scaled=ylim.scaled,
zlim.scaled=zlim.scaled, ...)
}
where pts is a dataframe with the three dimensional points (x,y,z) and the
conditioning variable g.
Mark
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] barchart (lattice) with text labels
Deepayan Sarkar wrote:
On 2/24/07, Mark and Heather Lyman [EMAIL PROTECTED] wrote:
I would like to place the value for each bar in barchart (lattice) at
the top of each bar. Something like the following code produces.
library(lattice)
mypanelfunc - function(x, y, ...)
{
panel.barchart(x, y, ...)
panel.text(x, y, labels=as.character(round(x,2)), ...)
}
myprepanelfunc - function(x, y, ...) list(xlim=c(0, max(x)+.1))
mydata - expand.grid(a=factor(1:5), b=factor(1:3), c=factor(1:2))
mydata$x - runif(nrow(mydata))
barchart(a~x|b, mydata, groups=c, panel=mypanelfunc,
prepanel=myprepanelfunc, adj=c(-0.1,0.5))
However, I cannot figure out how to shift the values to correspond with
their respective grouped bar.
You should look at panel.barchart and try to reproduce the
calculations done there.
Deepayan
This is the panel function that I ended up using. As Deepayan suggested,
I borrowed heavily from panel.barchart.
mypanelfunc - function(x, y, groups, box.ratio, ...)
{
panel.barchart(x, y, groups=groups, box.ratio=box.ratio, ...)
origin - current.panel.limits()$xlim[1]
groupSub - function(groups, subscripts, ...) groups[subscripts]
groups - as.numeric(groupSub(groups, ...))
vals - sort(unique(groups))
nvals - length(vals)
height - box.ratio/(1+ nvals * box.ratio)
for (i in unique(y)) {
ok - y == levels(y)[i]
nok - sum(ok, na.rm = TRUE)
panel.text(x = x[ok], y = (i + height * (groups[ok] - (nvals + 1)/2)),
labels=as.character(signif(x[ok],2)), ...)
}
}
Mark Lyman
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] sapply again return value
Antje niederlein-rstat at yahoo.de writes: Hello, I used an sapply to get some data back (s - sapply(...) ). The output of s would then deliver something like this: B06_lamp.csv C06_lamp.csv D06_lamp.csv [1,] NULL NULL Numeric,512 [2,] NULL NULL Numeric,512 [3,] NULL NULL 2 mode(s) [1] list dim(s) [1] 3 3 Now, I'd like to remove the columns which contain NULL (it's alway the whole column). How can I do this??? Antje As long as it is always the whole column, you can just test the first row, like this: s-matrix(list(NULL,NULL,NULL,NULL,NULL,NULL,numeric(512),numeric (512),2),ncol=3) s[,!apply(s,2,sapply,is.null)[1,],drop=FALSE] [,1] [1,] Numeric,512 [2,] Numeric,512 [3,] 2 If you would like to make sure that the whole column is NULL, you could do something like the following: s[,!apply(apply(s,2,sapply,is.null),2,sum)==nrow(s),drop=FALSE] [,1] [1,] Numeric,512 [2,] Numeric,512 [3,] 2 I use the drop=FALSE here for display purposes, but you may want to leave it off depending on desired format. Mark __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] LM Model
Simon P. Kempf simon.kempf at web.de writes: [SNIP] But I want to skip the lm function and specify my own regression equation RENT= 15 -0.15*AGE1 and then use the predict.lm function. However, in order to use the predict.lm function I need an object of class lm. Is there any way to do so? Or maybe somebody has another solution? Thanks in advance, Simon Here is one way. Take a look at help(offset), help(lm), and help(lm.predict). xx - runif(30) yy - rnorm(30) mydata-data.frame(xx,yy) lm(yy~offset(15*rep(1,30))+offset(-0.15*xx)-1,mydata) Mark Lyman __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Using variable names in for loops - Generating plot s semi-automatically from a series of variables Partly solved
Anthony Staines anthony.staines at gmail.com writes:
For each
disease, I want to pull out the column of data containing the word
'Male' and plot this against age, and then add a line to the plot for
the corresponding column containing 'Female'.
--
attach(data)
Diseases - c(Cardiovascular disease,Road Traffic Injury, ...
,All causes)
Male - names(data)[grep(Male,names(data))]
Female - names(data)[grep(Female,names(data))]
#Disease contains disease labels in the correct order, and Male and
Female now hold the (correct) variables.
for (i in seq(1,length(Diseases)))
{
jpeg(paste(Diseases[i],.jpg)) #This works fine!
plot(Male[i]~Age) #This does not
lines(Female[i]~Age)
dev.off()
}
detach(data)
[SNIP]
There are a few ways you can do this. Using plot.formula like you do here you
can use this:
mydata - as.data.frame(matrix(runif(300),30,10))
attach(mydata)
plot(formula(paste(mynames[1],mynames[2],sep='~')))
A way to do this without a formula would be:
plot(eval(parse(text=mynames[1])),eval(parse(text=mynames[2])))
detach(mydata)
Take a look at the help files for eval and parse. I still do not have a firm
grasp on how to use them and other related functions, like substitute, but
what I have been able figure out has been very useful.
Mark Lyman
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] plotting derived values not raw
James Root jcroot at gmail.com writes: I am trying to plot the mean and standard error of three separate conditions. For various reasons, I do not have access to the raw data from which the mean and error were derived and would like to make error bar plots utilizing only the actual mean and standard error values. Is there a way to do this in R? Thanks for any help in advance. james I believe that xYplot in the Hmisc package will do what you want Mark Lyman __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] applying cor.test in two different but matched dataframes
Gabriele Stocco stoccog at units.it writes: Dear R help, I have two dataframes with the same number of rows and columns. Each column is for a patient, each row for a variable. The columns and rows are matched in the two dataframes, so that each patient (column) and each variable (row) has the same position in the two dataframes. The values are different since each dataframe reports data from a different way to measure the same variables in the sampe patients. How can I make a cor.test, applying the function row by row in the two different dataframes, possibly without merging them? Thanks, Gabriele Stocco University of Trieste Try: # Create Data Frames With Rows as Variables x-as.data.frame(t(as.data.frame(matrix(rnorm(100),10 y-as.data.frame(t(as.data.frame(matrix(rnorm(100),10 # Convert to Lists So that Each Element of List is a variable t.x-as.list(as.data.frame(t(x))) t.y-as.list(as.data.frame(t(y))) # Use mapply to apply function to the two lists; see ?mapply mapply(cor.test,t.x,t.y) Mark Lyman __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Problem with fitted.value function
I don't know how to use the fitted values function with a given function and given input-variables but yet unknown result-values. Take a look at the predict function, ?predict. Mark Lyman __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] update index in
Kim Milferstedt milferst at uiuc.edu writes: Hello, I have a time series of data as a data.frame. Occasionally there is one or more days missing (e.g. data available for days 2, 3, 4, 8, 9, 10 -- missing days between 4 and 8). The experimental time information can be found in the 2nd column of data. I would like to have a continuous time line with one time point per day. Therefore I try to insert lines for the missing days that contain zeros for the data categories just to fill the columns. SNIP Thanks already, Kim SNIP I believe this will also do what you want: days-c(1:10)[-5:-7] xx-rnorm(7) data-data.frame(xx,days) new.data-merge(data,data.frame(days=1:10),all.y=TRUE) It usually is not a good idea to use zeroes as placeholders for missing values. Mark Lyman __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Plotmath expression
Christos Hatzis christos at nuverabio.com writes: Hello, I've been trying to plot a subscript in a text formula using plotmath but I haven't been able to do so. In my example below I would like the text label to show X[min] = 10.1 +/- 5.5 Here is the code: ll - c(x=10.1, sde=5.5) plot(1:10) text(x=9, y=2, pos=2, expression(paste(X[min], =, paste(ll, collapse=+/- How about the following? text(x=9, y=2, pos=2,substitute(X[min]==x%+-%sde,as.list(ll))) Mark Lyman __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] write data to pdf
Franco Mendolia franco.mendolia at gmx.de writes: Hello! Is there a possibility in R to save data in pdf-format? I do not want to save a plot but some lines of simple text. Regards, Franco Mendolia You could also use pdf() and textplot() in the gplots package Mark Lyman __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Incompatability of lme4 and Matrix packages
I just installed R 2.4 on a windows machine and installed the lme4 package version 0.9975-1. Unfortunately this version requires the 0.9975-2 version of the Matrix package which I cannot seem to find. It seems the newest version I can find is 0.9975-1. Is there somewhere I could find this newer version. Mark Lyman __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] How can I delete components in a column ?
Yen Ngo yen_h_ngo at yahoo.se writes: Hi all R-helpers, i am a new R-user and have problem with deleting some components in a column. I have a dataset like Name Idx empty 2 empty 3 anone2 bnone3 d none 2 ad cfh 4 bf cdt 5 empty 2 empty 2 gf cdh 4 d none 5 and want to eliminate all components that have id=none and empty . The remaining data should be Take a look at the subset function. Mark Lyman __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Help on apply() function
Anusha anu1978 at hotmail.com writes: Respected Sir/Madam, I have a problem with apply function. I have to two matrices of dimension of one column but n rows. I have to check whether one matrix is greater than other by going thru each row (ie) using if condition to check one matrix with another matrix. I like to use apply() function to this approach. That is apply function between two matrices. I searched for examples online but I couldn't find any. I don't know how to loop thru the matrices. Please help with this problem. Any help is appreciated. My mailid is monkponu at yahoo.com. Thanks, Anusha. You can use the functions all.equal with the function isTRUE (see ?all.equal) to check if two objects are nearly equal (within a certain tolerance). Or you can use identical (see ?identical) to check if they are exactly the same. See the examples in the help for identical. Mark Lyman __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] How to iteratively extract elements out of a list
lapply(m,function(x)x[x2]) [[1]] [1] 3 4 [[2]] [1] 4 5 [[3]] [1] 4 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Type II and III sum of square in Anova (R, car package)
1. First of all, more general question. Standard anova() function for lm() or aov() models in R implements Type I sum of squares (sequential), which is not well suited for unbalanced ANOVA. Therefore it is better to use Anova() function from car package, which was programmed by John Fox to use Type II and Type III sum of squares. Did I get the point? 2. Now more specific question. Type II sum of squares is not well suited for unbalanced ANOVA designs too (as stated in STATISTICA help), therefore the general rule of thumb is to use Anova() function using Type II SS only for balanced ANOVA and Anova() function using Type III SS for unbalanced ANOVA? Is this correct interpretation? 3. I have found a post from John Fox in which he wrote that Type III SS could be misleading in case someone use some contrasts. What is this about? Could you please advice, when it is appropriate to use Type II and when Type III SS? I do not use contrasts for comparisons, just general ANOVA with subsequent Tukey post-hoc comparisons. There are many threads on this list that discuss this issue. Not being a great statistician myself, I would suggest you read through some of these as a start. As I understand, the best philosophy with regards to types of sums of squares is to use the type that tests the hypothesis you want. They were developed as a convenience to test many of the hypotheses a person might want automatically, and put it into a nice, neat little table. However, with an interactive system like R it is usually even easier to test a full model vs. a reduced model. That is if I want to test the significance of an interaction, I would use anova(lm.fit2, lm.fit1) where lm.fit2 contains the interaction and lm.fit2 does not. The anova function will return the appropriate F-test. The danger with worrying about what type of sums of squares to use is that often we do not think about what hypotheses we are testing and if those hypotheses make sense in our situation. Mark Lyman __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Specify Z matrix with lmer function
Is there a way to specify a Z matrix using the lmer function, where the model is written as y = X*Beta + Z*u + e? I am trying to reproduce smoothing methods illustrated in the paper Smoothing with Mixed Model Software my Long Ngo and M.P. Wand. published in the /Journal of Statistical Software/ in 2004 using the lme4 and Matrix packages. The code and data sets used can be found at http://www.jstatsoft.org/v09/i01/. There original code did not work for me without slight modifications here is the code that I used with my modifications noted. x - fossil$age y - 10*fossil$strontium.ratio knots - seq(94,121,length=25) n - length(x) X - cbind(rep(1,n),x) Z - outer(x,knots,-) Z - Z*(Z0) # I had to create the groupedData object with one group to fit the model I wanted grp - rep(1,n) grp.dat-groupedData(y~Z|grp) fit - lme(y~-1+X,random=pdIdent(~-1+Z),data=grp.dat) I would like to know how I could fit this same model using the lmer function. Specifically can I specify a Z matrix in the same way as I do above in lme? Thanks, Mark __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Post Hoc Groupings
Looking at the errors your code produces, it looks like you need to make Dock and Slip factors. dock_2004_data$Dockf-factor(dock_2004_data$Dock) dock_2004_data$Slipf-factor(dock_2004_data$Slip) rich.aov - aov(X.open ~ Dockf*Slipf, data=dock_2004_data) TukeyHSD(rich.aov, c(Dockf, Slipf)) Jarrett Byrnes wrote: Indeed, the following works as well On Oct 26, 2005, at 5:23 PM, P Ehlers wrote: fm1 - aov(breaks ~ wool*tension, data = warpbreaks) TukeyHSD(fm1, c(wool,tension, wool:tension)) However, when working with my own dataset, I get the following errors. I have some inkling this may be due to a slightly unbalanced sample size, but am uncertain of this. rich.aov - aov(X.open ~ Dock*Slip, data=dock_2004_data) TukeyHSD(rich.aov, c(Dock, Slip)) Error in rep.int(n, length(means)) : unimplemented type 'NULL' in 'rep' In addition: Warning messages: 1: non-factors ignored: Slip in: replications(paste(~, xx), data = mf) 2: non-factors ignored: Dock, Slip in: replications(paste(~, xx), data = mf) I am pleased to know that these errors are not quite the __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Bug in lmer?
I am using version 0.98-7 of the Matrix package. I used the RGui
Install Packages... menu option to get the lme4 package from CRAN and
this version of the Matrix was automatically downloaded as well.
Martin Maechler wrote:
Mark == Mark Lyman [EMAIL PROTECTED]
on Thu, 29 Sep 2005 14:44:38 -0600 writes:
Mark I am relatively new to R so I am not confident enough in what I am
doing
Mark to be certain this is a bug.
Mark I am running R 2.1.1 on a Windows XP
Mark machine and the lme4 package version 0.98-1.
lme4 nowadays is heavily based on Matrix which version are you
using there?
Mark The following code fits the model I want using the
Mark nlme package version 3.1-60.
.. {see a script at the end}
Mark The problem is that when I try fitting the model using
Mark the lmer function with the following code:
Mark lmer(adg~trt+(1|loc)+(1|block:loc)+(1|loc:trt),mltloc)
Mark I get this message from Windows and R closes.
Mark R for Windows GUI front-end has encountered a problem and needs
to
Mark close. We are sorry for the inconvenience.
That definitely means there is a bug.
The question is *where* the bug is: lme4, Matrix, R, Windows.
One first thing coming to mind is a mismatch of lme4 and Matri
Mark This same code works on a Macintosh. So it doesn't
Mark seem that I have made an error in my code.
correct; I can also fit the model nice and quickly on Linux,
and summary() confirms the same fit {with the usual problem of
different estimates for the degrees of freedoms 'df'}.
So currently the bug only shows on the Windows platform.
Could it be that you have a mismatching package Matrix version
there, but not on the Mac?
Mark Also if anyone of the random effect terms is removed there is
Mark no problem. Is this something that is being looked at?
not yet, AFAIK.
Mark Or I have I made a mistake somewhere? I have included
Mark the data that I am using below.
I'm putting the data and an R script up for FTP,
so that you or others can run this ``from anywhere'' via
source(ftp://stat.ethz.ch/U/maechler/R/mltloc-ex.R;, echo = TRUE)
Maybe this helps diagnosis,
Martin Maechler
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] lmer random effect model matrix question
I have one fixed effect, sor, with two levels. I have eight lots and three wafers from each lot. I have included the data below. I would like to fit a mixed model that estimates a covariance parameter for wafer, which is nested in lot, and two covariance parameters for lot, one for each level of sor. The following command fits the model that I want, except for it estimates the correlation between the two covariance parameters for lot. Is there anyway to make R not estimate this correlation? Thank you. lmer(y~sor+(sor-1|lot)+(1|wafer:lot),wafer) For those familiar with proc mixed the following SAS code fits the model that I want: proc mixed scoring=4; class sor lot wafer site; model y= sor/ddfm=satterth; random lot(sor)/group=sor; random wafer(lot); run; sor lot wafer sitey 11 1 11 2006 21 1 12 1999 31 1 13 2007 41 1 21 1980 51 1 22 1988 61 1 23 1982 71 1 31 2000 81 1 32 1998 91 1 33 2007 10 1 2 11 1991 11 1 2 12 1990 12 1 2 13 1988 13 1 2 21 1987 14 1 2 22 1989 15 1 2 23 1988 16 1 2 31 1985 17 1 2 32 1983 18 1 2 33 1989 19 1 3 11 2000 20 1 3 12 2004 21 1 3 13 2004 22 1 3 21 2001 23 1 3 22 1996 24 1 3 23 2004 25 1 3 31 1999 26 1 3 32 2000 27 1 3 33 2002 28 1 4 11 1997 29 1 4 12 1994 30 1 4 13 1996 31 1 4 21 1996 32 1 4 22 2000 33 1 4 23 2002 34 1 4 31 1987 35 1 4 32 1990 36 1 4 33 1995 37 2 5 11 2013 38 2 5 12 2004 39 2 5 13 2009 40 2 5 21 2023 41 2 5 22 2018 42 2 5 23 2010 43 2 5 31 2020 44 2 5 32 2023 45 2 5 33 2015 46 2 6 11 2032 47 2 6 12 2036 48 2 6 13 2030 49 2 6 21 2018 50 2 6 22 2022 51 2 6 23 2026 52 2 6 31 2009 53 2 6 32 2010 54 2 6 33 2011 55 2 7 11 1984 56 2 7 12 1993 57 2 7 13 1993 58 2 7 21 1992 59 2 7 22 1992 60 2 7 23 1990 61 2 7 31 1996 62 2 7 32 1993 63 2 7 33 1987 64 2 8 11 1996 65 2 8 12 1989 66 2 8 13 1996 67 2 8 21 1997 68 2 8 22 1993 69 2 8 23 1996 70 2 8 31 1990 71 2 8 32 1989 72 2 8 33 1992 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Bug in lmer?
I am relatively new to R so I am not confident enough in what I am doing to be certain this is a bug. I am running R 2.1.1 on a Windows XP machine and the lme4 package version 0.98-1. The following code fits the model I want using the nlme package version 3.1-60. mltloc$loc - factor(mltloc$loc) mltloc$block - factor(mltloc$block) mltloc$trt - factor(mltloc$trt) Mltloc - groupedData(adg~trt|loc,mltloc) plot(Mltloc) plot(Mltloc,outer=~trt) lme(adg~trt, random=pdBlocked(list(pdIdent(~1),pdIdent(~block-1),pdIdent(~trt-1))),Mltloc) The problem is that when I try fitting the model using the lmer function with the following code: lmer(adg~trt+(1|loc)+(1|block:loc)+(1|loc:trt),mltloc) I get this message from Windows and R closes. R for Windows GUI front-end has encountered a problem and needs to close. We are sorry for the inconvenience. This same code works on a Macintosh. So it doesn't seem that I have made an error in my code. Also if anyone of the random effect terms is removed there is no problem. Is this something that is being looked at? Or I have I made a mistake somewhere? I have included the data that I am using below. X obs loc block trt adg fe 1 1 3 A 1 3 3.16454 7.1041 2 2 4 A 1 4 3.12500 6.6847 3 3 6 A 1 2 3.15944 6.8338 4 4 7 A 1 1 3.25000 6.5254 5 5 9 A 2 2 2.71301 8.2505 6 6 10 A 2 1 3.20281 7.5922 7 7 12 A 2 3 3.02423 7.3894 8 8 16 A 2 4 2.87245 7.4604 9 9 19 A 3 1 2.68878 8.2785 10 10 20 A 3 2 2.86862 7.9470 11 11 21 A 3 3 2.89923 7.9739 12 12 22 A 3 4 3.02806 8.4331 13 13 25 B 1 3 2.18131 6.6691 14 14 27 B 1 4 2.51914 5.6281 15 15 28 B 1 2 1.88739 7.0723 16 16 31 B 1 1 2.34685 6.0295 17 17 33 B 2 4 2.45608 5.6195 18 18 35 B 2 1 2.25225 6.3978 19 19 36 B 2 3 2.23649 6.1799 20 20 40 B 2 2 2.47523 5.9985 21 21 41 B 3 2 1.94200 7.2975 22 22 44 B 3 3 2.43243 6.4350 23 23 47 B 3 4 2.30180 6.3339 24 24 48 B 3 1 2.53378 6.1564 25 25 50 C 1 4 2.96014 7.5110 26 26 51 C 1 2 3.23551 7.4762 27 27 54 C 1 3 3.24638 7.2063 28 28 56 C 1 1 3.04710 7.6389 29 29 58 C 2 3 3.26449 7.5466 30 30 59 C 2 2 2.71377 9.0895 31 31 61 C 2 1 3.06522 7.8723 32 32 62 C 2 4 2.71739 8.2318 33 33 66 C 3 4 3.03623 7.9426 34 34 67 C 3 3 3.10507 8.4608 35 35 69 C 3 1 3.16304 8.5549 36 36 70 C 3 2 3.02899 8.5038 37 37 74 D 1 1 2.49164 9.5758 38 38 77 D 1 3 2.51833 9.5121 39 39 79 D 1 2 2.35631 10.3264 40 40 80 D 1 4 2.30331 9.7715 41 41 81 D 2 3 2.72688 9.5628 42 42 83 D 2 2 2.59512 9.9414 43 43 85 D 2 1 2.56516 9.3887 44 44 88 D 2 4 2.91523 8.3158 45 45 89 D 3 3 2.57943 10.4416 46 46 90 D 3 4 2.98159 8.7710 47 47 93 D 3 2 2.35370 11.0148 48 48 94 D 3 1 2.21953 11.2417 49 49 99 E 1 3 2.84158 8.7886 50 50 100 E 1 4 2.65264 8.6946 51 51 102 E 1 2 2.47112 9.7143 52 52 103 E 1 1 2.89769 9.2401 53 53 105 E 2 2 2.57343 9.5353 54 54 106 E 2 1 2.99752 8.7538 55 55 110 E 2 4 2.95380 8.8210 56 56 112 E 2 3 3.08663 8.9427 57 57 114 E 3 1 2.72525 9.4308 58 58 115 E 3 2 2.75825 9.7721 59 59 116 E 3 3 3.08333 8.9010 60 60 117 E 3 4 3.12129 8.4852 61 61 122 F 1 1 3.20600 6.3983 62 62 123 F 1 2 2.89500 6.6569 63 63 125 F 1 4 3.36900 6.0821 64 64 126 F 1 3 3.12000 6.5349 65 65 130 F 2 2 3.19300 6.6729 66 66 131 F 2 1 3.29800 6.5488 67 67 133 F 2 4 3.09700 6.6598 68 68 135 F 2 3 3.38500 6.2998 69 69 139 F 3 3 3.44900 6.2849 70 70 140 F 3 2 3.05000 6.9957 71 71 141 F 3 1 3.43500 6.7302 72 72 143 F 3 4 3.60600 6.3827 73 73 145 G 1 2 2.58669 8.1394 74 74 147 G 1 1 3.17892 7.0972 75 75 148 G 1 4 2.95284 7.3140 76 76 151 G 1 3 3.17924 6.9430 77 77 154 G 2 4 2.62344 7.5150 78 78 155 G 2 3 2.64286 8.0237 79 79 157 G 2 1 3.12760 7.3169 80 80 160 G 2 2 2.54993 8.1957 81 81 163 G 3 4 2.58322 7.9687 82 82 164 G 3 3 2.84813 7.9284 83 83 166 G 3 2 2.69279 8.5303 84 84 167 G 3 1 3.14424 7.3564 85 85 169 H 1 3 3.39974 6.5945 86 86 173 H 1 2 3.12370 6.7530 87 87 175 H 1 4 3.17969 6.4279 88 88 176 H 1 1 3.70052 6.4830 89 89 177 H 2 2 2.95192 7.3809 90 90 179 H 2 3 3.44661 6.7929 91 91 182 H
