[R] Functions in lists or arrays?
I have a problem where I need to have a "driver" style R program that will
extend itself, given more
'source("extra.R")' style lines. IE: easy to modify by other people. The
problem becomes when I
would like to create an array or list of functions. Each function, possibly
hundreds of them, are really
created by other programs, generating an *.R file. So for example, if I
try:
> t <- list()
> t[1] <- function(b) { b*2 }
Error in t[1] <- function(b) { :
cannot coerce type 'closure' to vector of type 'list'
Similar errors for arrays, and anything else I can think of. I'm an R
neophite, which likely shows. The
only way I seem to be able to do the above, is to generate unique names for
each function, and add
each name into a list, sort of like this:
# Register ourselves
models <- cbind(models, "Nasa_PART_Rules")
bounds <- cbind(bounds, "Nasa_PART_Bounds")
Nasa_PART_Rules <- NULL
Nasa_PART_Bounds <- NULL
# Rules section
Nasa_PART_Rules <- rbind(Nasa_PART_Rules, c("Nasa_PART_R1", "F"))
Nasa_PART_R1 <- function(f) {
f[,"CYCLOMATIC_COMPLEXITY"] > 8 &
f[,"CYCLOMATIC_COMPLEXITY"] <= 60 &
f[,"LOC_TOTAL"] > 73
}
Nasa_PART_Bounds <- rbind(Nasa_PART_Bounds, c("Nasa_PART_B1"))
Nasa_PART_B1 <- function(b) {
b["CYCLOMATIC_COMPLEXITY",0] <- 8
b["CYCLOMATIC_COMPLEXITY",1] <- 60
b["LOC_TOTAL",0] <- 73
}
#...
And then using something like this function:
# Dispatch a function from its name
dispatch <- function(f, x) {
eval(call(f, x))
}
to evaluate each rule over all the data rows:
# Read training+validation data
dat <- read.csv("jm1_data(training+validation).csv")
mat <- NULL
clt <- NULL
# Evaluate each rule against the dataset
for (i in models) {
# Get the rules for the model
rules <- eval(substitute(foo[,1], list(foo=as.name(i
cls <- eval(substitute(foo[,2], list(foo=as.name(i
res <- lapply(rules, dispatch, dat)
#...
Now, this seems way too uggly to me. Can someone give me a hand and/or
point me into a more sane
direction to explore?
One option I have thought of, is to get rid of the *_B?() functions and just
fill in a 3 dimensional array using
something like:
x <- NULL
dimnames(x) <- c(colnames(mat),colnames(dat), c("lbound","ubound"))
...
x["RULE_NAME_1", "DATA_COL_NAME_1", "lbound"] <- ...
...
But I'm not exactly sure how I would construct and/or add onto a global
array/etc extra dimnames, as I source
each generated *.R file.
Anyways, Not sure if I'm making much sense... thanks for any help,
-Toby.
[[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] Generalized 2D list/array/whatever?
I'm trying to figure out how I can get a generalized 2D list/array/matrix/whatever working. Seems I can't figure out how to make the variables the right type. I always seem to get some sort of error... out of bounds, wrong type, wrong dim, etc. Very confused... :) x[["some label", "some other index"]] <- 3 x[["some other label", "something else"]] <- 4 I don't know the indexes/label ahead of time... they get generated... Any thoughts? -Toby. [[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] Generalized 2D list/array/whatever?
I must be going about my idea really backwards. I have functions like so:
Nasa_PART_Bounds <- rbind(Nasa_PART_Bounds, c("Nasa_PART_B1"))
Nasa_PART_B1 <- function(l,u) {
l[["Nasa_PART_B1", "CYCLOMATIC_COMPLEXITY"]] <- 8
u[["Nasa_PART_B1", "CYCLOMATIC_COMPLEXITY"]] <- 60
l[["Nasa_PART_B1", "LOC_TOTAL"]] <- 73
}
Where I add the function names of each "bounds" function to a list,
which later on get called with two variables.
dispatch2 <- function(f, x, y) {
eval(call(f, x, y))
}
for (i in bounds) {
bb <- eval(substitute(foo[,1], list(foo=as.name(i
rr <- lapply(bb, dispatch2, lb, ub)
}
However, I don't know how to do this in a better "R" way. I think my C and
pearl experience is showing through... :(
-Toby.
On Thu, Apr 23, 2009 at 8:20 PM, Gabor Grothendieck wrote:
> Matrix made from a list:
>
> m <- list(sin, 1:3, letters[1:3], expression(a+b))
> dim(m) <- c(2, 2)
> dimnames(m) <- list(letters[1:2], LETTERS[1:2])
> class(m) # matrix
>
> or
>
> M <- structure(list(sin, 1:3, letters[1:3], expression(a+b)), .Dim = c(2,
> 2),
> .Dimnames = list(c("a", "b"), c("A", "B")))
> class(M) # matrix
>
> On Thu, Apr 23, 2009 at 10:03 PM, Toby
> wrote:
> > I'm trying to figure out how I can get a generalized 2D
> > list/array/matrix/whatever
> > working. Seems I can't figure out how to make the variables the right
> > type. I
> > always seem to get some sort of error... out of bounds, wrong type, wrong
> > dim, etc.
> > Very confused... :)
> >
> > x[["some label", "some other index"]] <- 3
> > x[["some other label", "something else"]] <- 4
> >
> > I don't know the indexes/label ahead of time... they get generated...
> Any
> > thoughts?
> >
> > -Toby.
> >
> >[[alternative HTML version deleted]]
> >
> > __
> > R-help@r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> > and provide commented, minimal, self-contained, reproducible code.
> >
>
[[alternative HTML version deleted]]
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
[R] function for between-groups covariance matrix?
Hi Is anyone aware if there is a function already available that calculates the between-groups COvariance matrix, say in a discriminant analysis setting, or in a manova setting. Maybe as a helper function to some other major function. Otherwise I would have to program it myself (probably resulting in some awkward code performing horribly on this large dataset). Thanks Toby __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] nls() fails on a simple exponential fit, when lm() gets it right?
Dear R-help,
Here's a simple example of nonlinear curve fitting where nls seems to get
the answer wrong on a very simple exponential fit (my R version 2.7.2).
Look at this code below for a very basic curve fit using nls to fit to (a)
a logarithmic and (b) an exponential curve. I did the fits using
self-start functions and I compared the results with a more simple fit
using a straight lm() command.
The logarithmic fit works 100% correctly below. The problem is with the
exponential fit: the nls command gives the wrong values and I have
completely failed to see why this should happen. I can't see any misake in
my self-start function, so why the mismatch?? Even worse, if I give nls a
trivial fit by removing the "#&#&" on line 6, it fails to converge at all
when the 'simpler' method with lm() easily gives the right answer (r=0.03
and A=5).
I did the same exp and ln fits using MS Excel and the numbers returned are
exactly as for the 'simpler ways', so nls is definitely in the wrong, but
the self-start is almost identical to the one for the logarithmic fit,
which works perfectly
I can't see my mistake at all. Can anyone help??
Thanks,
Toby Marthews
traceson=FALSE;maxint=1;minstepfactor=2^(-30);tolerance=1e-7
xtxt="DBH in mm";ytxt="Height in m"
dbh=c(0.9,1.0,1.1,1.2,4.8,4.9,5.0,5.1,8.9,9.0,9.1,9.2,11.8,11.9,12.0,12.1,14.9,15.1,15.2,15.5)
height=c(5.770089,5.154794,4.47,6.356770,14.849109,14.973146,15.423816,14.865038,21.335985,20.477624,20.915823,21.886004,23.714724,24.182210,23.954929,23.784659,25.334360,25.411320,26.218614,25.612478)
#&#&#height=5*exp(0.03*dbh)
logarithmicInit=function(mCall,LHS,data) {
xy=sortedXyData(mCall[["x"]],LHS,data)
aaa=0.01 #Just guess aaa=0.01 to start the fit
bbb=min(xy$y) #Use the minimum y value as an initial guess for bbb
value=c(aaa,bbb)
names(value)=mCall[c("aaa","bbb")]
return(value)
}
fmla=as.formula("~(aaa*log(x)+bbb)")
SSlogarithmic=selfStart(fmla,initial=logarithmicInit,parameters=c("aaa","bbb"))
exponenInit=function(mCall,LHS,data) {
xy=sortedXyData(mCall[["x"]],LHS,data)
r=0.01 #Just guess r=0.01 to start the fit
A=min(xy$y)#Use the minimum y value as an initial guess for A
value=c(r,A)
names(value)=mCall[c("r","A")]
return(value)
}
fmla=as.formula("~(A*exp(r*x))")
SSexponen=selfStart(fmla,initial=exponenInit,parameters=c("r","A"))
cat("Logarithmic fit (here, the self-start and the 'simpler' way match):\n")
tmp=getInitial(height~SSlogarithmic(dbh,aaa,bbb),data=data.frame(dbh,height))
cat("(Starting values for the logarithmic fit were: aaa=",tmp[1],"and
bbb=",tmp[2],")\n")
modelfit=nls(height~SSlogarithmic(x=dbh,aaa,bbb),trace=traceson,control=list(maxiter=maxint,tol=tolerance,minFactor=minstepfactor))
bfaaa=summary(modelfit)$coefficients[1];bfbbb=summary(modelfit)$coefficients[2]
cat("Best logarithmic fit is (",ytxt,")=(aaa*log(",xtxt,")+bbb) with
parameters aaa=",bfaaa,"and bbb=",bfbbb,"\n")
#Produces: Best logarithmic fit is ( Height in m )=(aaa*log( DBH in mm
)+bbb) with parameters aaa= 7.551666 and bbb= 4.594367
lnfit=lm(height~log(dbh)) #This is doing the ln fit without a self-start
function
cat("Done another, simpler, way, the best logarithmic fit has parameters:
aaa=",lnfit$coefficients[2],"and bbb=",lnfit$coefficients[1],"\n")
#Produces: Done another, simpler, way, the best logarithmic fit has
parameters: aaa= 7.551666 and bbb= 4.594367
cat("Exponential fit (here, the self-start and the 'simpler' way DON'T
match):\n")
tmp=getInitial(height~SSexponen(dbh,r,A),data=data.frame(dbh,height))
cat("(Starting values for the exponential fit: r=",tmp[1],"and
A=",tmp[2],")\n")
modelfit=nls(height~SSexponen(x=dbh,r,A),trace=traceson,control=list(maxiter=maxint,tol=tolerance,minFactor=minstepfactor))
bfr=summary(modelfit)$coefficients[1];bfA=summary(modelfit)$coefficients[2]
cat("Best exponential fit is (",ytxt,")=(A*exp(r*(",xtxt,"))) with
parameters r=",bfr,"and A=",bfA,"\n")
#Produces: Best exponential fit is ( Height in m )=(A*exp(r*( DBH in mm
))) with parameters r= 0.07134055 and A= 9.47907
expfit=lm(log(height)~dbh) #This is doing the exp fit without a self-start
function
cat("Done another, simpler, way, the best exponential fit has parameters:
r=",expfit$coefficients[2],"and A=",exp(expfit$coefficients[1]),"\n")
#Produces: Done another, simpler, way, the best exponential fit has
parameters: r= 0.1028805 and A= 6.75045
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
[R] extracting values conditonal on other values
Dear R helpers,
I have a dataframe (test1) containing the time of sunrise and sunset for each
day of the year
for 3 years. I have another dataframe (test2) containing measurements that are
taken every
15 minutes, 24/7. I would like to extract all rows from test2 that occur
between sunrise and
sunset for the appropriate date. Can you suggest a good vectorized way to do
this? Keep in
mind that the sunrise/sunset dataframe has 1 row for each day, and the
measurement
dataset has 96 rows for each day. I'm hoping not to have to match strings...
The times (test1$rise, test1$set, and test2$time) in the example are rather
ugly since I wasn't
sure how to generate a random hourly time series. I would also use a standard
date format
for the real thing but, again, wasn't sure how to generate dates for this
example.
Example data:
test1 <- data.frame(year = gl(3, 30, 90, labels = c("2006", "2007", "2008")),
month = gl(3, 10, 30, labels = c("1", "2", "3")),
day = rep(c(21:30, 19:28, 22:31),3),
rise = as.integer(runif(90, 700, 750)),
set = as.integer(runif(90, 1630,1745)))
test2 <- data.frame(year = gl(3, 2880, 8640, labels = c("2006", "2007",
"2008")),
month = gl(3, 96, 288, labels = c("1", "2", "3")),
day = rep(c(21:30, 19:28, 22:31),3, each = 96),
time = 100*rep(seq(, 23.75,by= .25),90),
temp = runif(8640, -5, 15))
Thank you in advance,
Toby
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Reshape
Hello, helpeRs,
I've been trying, unsuccessfully, to change a dataframe from long to
wide format using reshape (the original). I would appreciate it if
someone could demonstrate the correct syntax. The script below will
create a toy example. The new wide data should have a column name
for each unique entry in the "fame" column. Under each column
should be either the appropriate weight or na, if there is no match.
Thus, the end product should have column names:
valley plot trt 18w 16iso 12:0, etc.
Here is the toy script for the starting data in long form:
dat <- data.frame(fame = gl(4,1,10, labels = c( "18w", "16iso",
"12:0", "16w")), valley = gl(2,6,10, labels = c("H", "M")), plot =
gl(5, 2, 10), trt = gl(2,1,10, labels = c("e", "g")), weight =
1000*runif(10, 0, 1))
Thank you for your assistance.
Toby
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Reshape
Hello, helpeRs,
I've been trying, unsuccessfully, to change a dataframe from long to
wide format using reshape (the original). I would appreciate it if
someone could demonstrate the correct syntax. The script below will
create a toy example. The new wide data should have a column name
for each unique entry in the "fame" column. Under each column
should be either the appropriate weight or na, if there is no match.
Thus, the end product should have column names:
valley plot trt 18w 16iso 12:0, etc.
Here is the toy script for the starting data in long form:
dat <- data.frame(fame = gl(4,1,10, labels = c( "18w", "16iso",
"12:0", "16w")), valley = gl(2,6,10, labels = c("H", "M")), plot =
gl(5, 2, 10), trt = gl(2,1,10, labels = c("e", "g")), weight =
1000*runif(10, 0, 1))
Thank you for your assistance.
Toby
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Reshape
That does work, thank you. I didn't understand that the "fame"
column would be the time varying column.
Toby
On 28 Sep 2010 at 12:47, Michael Bedward wrote:
> Hi Toby,
>
> I think this should work...
>
> reshape(dat, v.names=c("weight"), idvar=c("valley", "plot", "trt"),
> timevar="fame", direction="wide")
>
> Michael
>
>
> On 28 September 2010 12:17, Toby Gass
> wrote:
> > Hello, helpeRs,
> >
> > I've been trying, unsuccessfully, to change a dataframe from long to
> > wide format using reshape (the original). I would appreciate it if
> > someone could demonstrate the correct syntax. The script below will
> > create a toy example. The new wide data should have a column name
> > for each unique entry in the "fame" column. Under each column
> > should be either the appropriate weight or na, if there is no match.
> > Thus, the end product should have column names:
> >
> > valley plot trt 18w 16iso 12:0, etc.
> >
> > Here is the toy script for the starting data in long form:
> >
> > dat <- data.frame(fame = gl(4,1,10, labels = c( "18w", "16iso",
> > "12:0", "16w")), valley = gl(2,6,10, labels = c("H", "M")), plot =
> > gl(5, 2, 10), trt = gl(2,1,10, labels = c("e", "g")), weight =
> > 1000*runif(10, 0, 1))
> >
> > Thank you for your assistance.
> >
> > Toby
> >
> > __
> > R-help@r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/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] POSIX turns into factor
Dear helpeRs, I am working with a dataframe that includes a column, "calendar", used for plotting time series. > class(dat$calendar) [1] "POSIXt" "POSIXlt" When I finish working, I save my data as a .csv file. When I read the file in again, "calendar" is always a factor > class(dat$calendar) [1] "factor" and I have to turn it back into a date-time object. Is this unavoidable when going back and forth from a .csv, or can I do something differently to retain the class? R 2.11.1 on Windows OS. Thank you. Toby __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Using anova() with glmmPQL()
Dear R HELP,
ABOUT glmmPQL and the anova command. Here is an example of a repeated-measures
ANOVA focussing on the way starling masses vary according to (i) roost
situation and (ii) time (two time points only).
library(nlme);library(MASS)
stmass=c(78,88,87,88,83,82,81,80,80,89,78,78,85,81,78,81,81,82,76,74,79,73,79,75,77,78,80,78,83,84,77,68,75,70,74,84,80,75,76,75,85,88,86,95,100,87,98,86,89,94,84,88,91,96,86,87,93,87,94,96,91,90,87,84,86,88,92,96,83,85,90,87,85,81,84,86,82,80,90,77)
roostsitu=factor(c("tree","tree","tree","tree","tree","tree","tree","tree","tree","tree","nest-box","nest-box","nest-box","nest-box","nest-box","nest-box","nest-box","nest-box","nest-box","nest-box","inside","inside","inside","inside","inside","inside","inside","inside","inside","inside","other","other","other","other","other","other","other","other","other","other","tree","tree","tree","tree","tree","tree","tree","tree","tree","tree","nest-box","nest-box","nest-box","nest-box","nest-box","nest-box","nest-box","nest-box","nest-box","nest-box","inside","inside","inside","inside","inside","inside","inside","inside","inside","inside","other","other","other","other","other","other","other","other","other","other"),levels=c("tree","nest-box","inside","other"))
mnth=factor(c(rep("Nov",times=40),rep("Jan",times=40)),levels=c("Nov","Jan"))
subjectnum=c(1:10,1:10,1:10,1:10,1:10,1:10,1:10,1:10)
subject=factor(paste(roostsitu,subjectnum,sep=""))
dataf=data.frame(mnth,roostsitu,subjectnum,subject,stmass)
lmeres=lme(fixed=stmass~mnth*roostsitu,random=~1|subject/mnth,na.action=na.exclude)
anova(object=lmeres,test="Chisq")
numDF denDF F-value p-value
(Intercept)136 31143.552 <.0001
mnth 13695.458 <.0001
roostsitu 33610.614 <.0001
mnth:roostsitu 336 0.657 0.5838
I can conclude from this that variation with both roost situation and month are
significant, but with no interaction term. So far so good. However, say I were
interested only in whether or not those starlings were heavier or lighter than
78g: seemingly, I could change my response variable and analyse like this -
stmassheavy=ifelse(stmass>78,1,0)
lmeres1=lme(fixed=stmassheavy~mnth*roostsitu,random=~1|subject/mnth,na.action=na.exclude,family=binomial)
anova(object=lmeres1,test="Chisq")
but I get errors doing that. After a certain amount of web searching, I find
that I'm supposed to use glmmPQL for this so I tried:
lmeres2=glmmPQL(fixed=stmassheavy~mnth*roostsitu,random=~1|subject/mnth,na.action=na.exclude,family=binomial)
anova(object=lmeres2,test="Chisq")
The glmmPQL command runs, but I get "Error in anova.glmmPQL(object = lmeres,
test = "Chisq") : 'anova' is not available for PQL fits". Looking into this, I
find that I am not supposed to use the anova command in conjunction with
glmmPQL (several posts from Brian Ripley
http://r.789695.n4.nabble.com/R-glmmPQL-in-2-3-1-td808574.html and
http://www.biostat.wustl.edu/archives/html/s-news/2002-06/msg00055.html and
http://www.mail-archive.com/r-help@stat.math.ethz.ch/msg46894.html ) even
though it appears that the earlier versions of glmmPQL did allow the anova
command to work (before ~2004).
However, I couldn't find any other way to run a repeated-measures ANOVA with
famiy=binomial. After a while longer on Google, I found a 'workaround' from
Spencer Graves (on
http://markmail.org/message/jddj6aq66wdidrog#query:how%20to%20use%20anova%20with%20glmmPQL+page:1+mid:jddj6aq66wdidrog+state:results
):
class(lmeres2)="lme"
anova(object=lmeres2,test="Chisq")
numDF denDF F-value p-value
(Intercept)136 182.84356 <.0001
mnth 136 164.57288 <.0001
roostsitu 336 17.79263 <.0001
mnth:roostsitu 336 3.26912 0.0322
Which does give me a result and tells me that the interaction term is
significant here. HOWEVER, on that link Douglas Bates told Spencer Graves that
this wasn't an approprate method.
I haven't found any other workarounds for this except some general advice that
I should move onto using the lmer command (which I can't do because I need to
get p-values for my fits and according to
https://stat.ethz.ch/pipermail/r-help/2006-May/094765.html I won't get those
from lmer).
My questions are: (1) Is lmer the only way to do a binomial repeated-measures
ANOVA in R? (which means that there's no way to do such an ANOVA in R without
losing the p-values) and (2) if I am supposed to be using glmmPQL for this
simple situation, what am I doing wrong?
Thanks very much for any help anyone can give me.
best,
Toby Marthews
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Statistical formulas
Hi Ilya, If you're looking for general information about statistical tests, etc., you'll probably need to buy yourself a textbook. There are online pages (http://homes.msi.ucsb.edu/~byrnes/rtutorial.html is a good one), but a good textbook is probably better. It's a bit old, but I still recommend Fowler et al. (1998) for those who are starting out with statistical testing: Fowler J, Cohen L & Jarvis P (1998). Practical statistics for field biology (2nd ed.). Wiley, Chichester, UK. Toby Marthews From: r-help-boun...@r-project.org [r-help-boun...@r-project.org] On Behalf Of ilya [flya...@gmail.com] Sent: 18 January 2011 16:23 To: r-help@r-project.org Subject: [R] Statistical formulas Hello! I hope, i'm doing right, writing here. I have a question about the R statistical formulas: were can I look, what are they? E. g. I'd like to know which exactly mathematical formula is used in Wilcoxon test, and others. Is there somewhere on the Internet information about that? Thanks a lot! Ilya. PS Sorry for any mistakes in my English, it is not my mother-tongue. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Using anova() with glmmPQL()
Hi Ben Bolker, Just to say thank you VERY much for the reply below and for taking the time to go through my code. I think you're absolutely right and I have been using the wrong formula completely. I have been trying out several more examples (am still struggling with this) and will submit any further posts to r-sig-mixed-models. Best, Toby From: r-help-boun...@r-project.org [r-help-boun...@r-project.org] On Behalf Of Ben Bolker [bbol...@gmail.com] Sent: 18 January 2011 16:15 To: r-h...@stat.math.ethz.ch Subject: Re: [R] Using anova() with glmmPQL() Toby Marthews ouce.ox.ac.uk> writes: > > Dear R HELP, > > ABOUT glmmPQL and the anova command. Here is an example of a > repeated-measures ANOVA focussing on the way > starling masses vary according to (i) roost situation and > (ii) time (two time points only). [snip] > lmeres=lme(fixed=stmass~mnth*roostsitu, > random=~1|subject/mnth,na.action=na.exclude) > anova(object=lmeres,test="Chisq") By the way, 'test="Chisq"' has no effect -- it's silently ignored -- in an anova run on an lme() object the F test is always used (as your output suggests). > numDF denDF F-value p-value > (Intercept)136 31143.552 <.0001 > mnth 13695.458 <.0001 > roostsitu 33610.614 <.0001 > mnth:roostsitu 336 0.657 0.5838 By the way, I'm a bit puzzled/suspicious about your setup here. You seem to measure each subject once in each month, so your grouping variable subject/month is confounded with the residual error. I would suggest random=~1|subject (which gives a very nearly identical ANOVA table but makes more sense). In this case, because you have only two months, you could do almost the same analysis by taking differences or averages of subjects across months: a paired t-test on subjects for the effect of month, a one-way ANOVA on subject averages for the effect of roost, and a one-way ANOVA on subject differences for the month x roost interaction. > I can conclude from this that variation with > both roost situation and month are significant, but with no > interaction term. So far so good. However, say I > were interested only in whether or not those starlings > were heavier or lighter than 78g: seemingly, I could > change my response variable and analyse like this - [snip lme with family="binomial"] > lmeres2=glmmPQL(fixed=stmassheavy~mnth*roostsitu, > random=~1|subject/mnth,na.action=na.exclude,family=binomial) > anova(object=lmeres2,test="Chisq") > > The glmmPQL command runs, but I get > "Error in anova.glmmPQL(object = lmeres, test = "Chisq") : > 'anova' is > not available for PQL fits". [snip] >However, I couldn't find any other way to run a repeated-measures ANOVA with famiy=binomial. After a while > longer on Google, I found a 'workaround' from Spencer Graves (on > http://markmail.org/message/jddj6aq66wdidrog#query:how%20to%20use%20anova%20with%20glmmPQL+page:1+mid:jddj6aq66wdidrog+state:results ): > > class(lmeres2)="lme" > anova(object=lmeres2,test="Chisq") > >numDF denDF F-value p-value > (Intercept)136 182.84356 <.0001 > mnth 136 164.57288 <.0001 > roostsitu 336 17.79263 <.0001 > mnth:roostsitu 336 3.26912 0.0322 > > Which does give me a result and tells me that the interaction term is significant here. HOWEVER, on that link > Douglas Bates told Spencer Graves that this wasn't an approprate method. > > I haven't found any other workarounds for this except some > general advice that I should move onto using the > lmer command (which I can't do because I need to get p-values > for my fits and according to > https://stat.ethz.ch/pipermail/r-help/2006-May/094765.html > I won't get those from lmer). > > My questions are: (1) Is lmer the only way to do > a binomial repeated-measures ANOVA in R? (which means that > there's no way to do such an ANOVA in R without losing the p-values) > and (2) if I am supposed to be using > glmmPQL for this simple situation, what am I doing wrong? > Hmm. This feels harder than it ought to be, and I've already spent longer on it than I should, so I'm going to give you some thoughts and encourage you to send this on to R-sig-mixed-models. I think the *best* way to do this would probably be to do the analogue of the paired/grouped tests I suggested above: for the month effect, test whether more individuals go from light to heavy or heavy to light, and so forth. Treating this as a binomial problem: you can use glmmPQL and use the Wald test (see wald.test in
[R] conditional selection of dataframe rows
Dear helpeRs, I have a dataframe (14947 x 27) containing measurements collected every 5 seconds at several different sampling locations. If one measurement at a given location is less than zero on a given day, I would like to delete all measurements from that location on that day. Here is a toy example: toy <- data.frame(CH = rep(3:5,3), DAY = c(rep(4,5), rep(5,4)), SLOPE = c(seq(0.2,0.6, .1),seq(0.2, -0.1, -0.1))) In this example, row 9 has a negative measurement for Chamber 5, so I would like to delete row 6, which is the same Chamber on the same day, but not row 3, which is the same chamber on a different day. In the full dataframe, there are, of course, many more days. Is there a handy R way to do this? Thank you for the assistance. Toby __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 selection of dataframe rows
Thank you all for the quick responses. So far as I've checked, Marc's solution works perfectly and is quite speedy. I'm still trying to figure out what it is doing. :) Henrique's solution seems to need some columns somewhere. David's solution does not find all the other measurements, possibly with positive values, taken on the same day. Thank you again for your efforts. Toby On 12 Aug 2010 at 14:32, Marc Schwartz wrote: > On Aug 12, 2010, at 2:24 PM, Marc Schwartz wrote: > > > On Aug 12, 2010, at 2:11 PM, Toby Gass wrote: > > > >> Dear helpeRs, > >> > >> I have a dataframe (14947 x 27) containing measurements collected > >> every 5 seconds at several different sampling locations. If one > >> measurement at a given location is less than zero on a given day, I > >> would like to delete all measurements from that location on that day. > >> > >> Here is a toy example: > >> > >> toy <- data.frame(CH = rep(3:5,3), DAY = c(rep(4,5), rep(5,4)), > >> SLOPE = c(seq(0.2,0.6, .1),seq(0.2, -0.1, -0.1))) > >> > >> In this example, row 9 has a negative measurement for Chamber 5, so I > >> would like to delete row 6, which is the same Chamber on the same > >> day, but not row 3, which is the same chamber on a different day. In > >> the full dataframe, there are, of course, many more days. > >> > >> Is there a handy R way to do this? > >> > >> Thank you for the assistance. > >> > >> Toby > > > > > > > > Not fully tested, but here is one possibility: > > > >> toy > > CH DAY SLOPE > > 1 3 4 0.2 > > 2 4 4 0.3 > > 3 5 4 0.4 > > 4 3 4 0.5 > > 5 4 4 0.6 > > 6 5 5 0.2 > > 7 3 5 0.1 > > 8 4 5 0.0 > > 9 5 5 -0.1 > > > > > >> subset(toy, ave(SLOPE, CH, DAY, FUN = function(x) any(x < 0)) == 0) > > CH DAY SLOPE > > 1 3 4 0.2 > > 2 4 4 0.3 > > 3 5 4 0.4 > > 4 3 4 0.5 > > 5 4 4 0.6 > > 7 3 5 0.1 > > 8 4 5 0.0 > > > This can actually be slightly shortened to: > > > subset(toy, !ave(SLOPE, CH, DAY, FUN = function(x) any(x < 0))) > CH DAY SLOPE > 1 3 4 0.2 > 2 4 4 0.3 > 3 5 4 0.4 > 4 3 4 0.5 > 5 4 4 0.6 > 7 3 5 0.1 > 8 4 5 0.0 > > > HTH, > > Marc > __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 selection of dataframe rows
Hi, I do want to look only at slope. If there is one negative slope measurement for a given day and a given chamber, I would like to remove all other slope measurements for that day and that chamber, even if they are positive. On one day, I will have 20 slope measurements for each chamber. If one is negative, I would like to delete the other 19 for that chamber on that day, even if they are positive. I have measurements for every day of the year, for 4 years and multiple chambers. I know I could make some awful nested loop with a vector of day and chamber numbers for each occurrence of a negative slope and then run that against the whole data set but I hope not to have to do that. Here is the rationale, if that helps. These are unattended outdoor chambers that measure soil carbon efflux. When the numbers go negative during part of the day but otherwise look normal, it usually means a plant has sprouted in the chamber and is using the carbon dioxide. That means the measurements are all lower than they should be and I need to discard all measurements collected on that day, whether positive or negative. It might have been a little clearer if I'd make the toy dataframe a bit larger. Thanks again for the assistance. Toby On 12 Aug 2010 at 16:39, David Winsemius wrote: > > On Aug 12, 2010, at 4:06 PM, Toby Gass wrote: > > > Thank you all for the quick responses. So far as I've checked, > > Marc's solution works perfectly and is quite speedy. I'm still > > trying to figure out what it is doing. :) > > > > Henrique's solution seems to need some columns somewhere. David's > > solution does not find all the other measurements, possibly with > > positive values, taken on the same day. > > I assumed you only wanted to look at what appeared to be a data > column, SLOPE. If you want to look at all columns for negatives then > try: > > toy[ which( apply(toy, 1, function(x) all(x >= 0)) ), ] # or > toy[ apply(toy, 1, function(x) all(x >= 0)) , ] > > This is how they differ w,r,t, their handling of NA's. > > > toy[3,2] <- NA > > toy[ apply(toy, 1, function(x) all(x >= 0)) , ] > CH DAY SLOPE > 1 3 4 0.2 > 2 4 4 0.3 > NA NA NANA > 4 3 4 0.5 > 5 4 4 0.6 > 6 5 5 0.2 > 7 3 5 0.1 > 8 4 5 0.0 > > toy[ which(apply(toy, 1, function(x) all(x >= 0)) ), ] >CH DAY SLOPE > 1 3 4 0.2 > 2 4 4 0.3 > 4 3 4 0.5 > 5 4 4 0.6 > 6 5 5 0.2 > 7 3 5 0.1 > 8 4 5 0.0 > > > > > > Thank you again for your efforts. > > > > Toby > > > > On 12 Aug 2010 at 14:32, Marc Schwartz wrote: > > > >> On Aug 12, 2010, at 2:24 PM, Marc Schwartz wrote: > >> > >>> On Aug 12, 2010, at 2:11 PM, Toby Gass wrote: > >>> > >>>> Dear helpeRs, > >>>> > >>>> I have a dataframe (14947 x 27) containing measurements collected > >>>> every 5 seconds at several different sampling locations. If one > >>>> measurement at a given location is less than zero on a given day, I > >>>> would like to delete all measurements from that location on that > >>>> day. > >>>> > >>>> Here is a toy example: > >>>> > >>>> toy <- data.frame(CH = rep(3:5,3), DAY = c(rep(4,5), rep(5,4)), > >>>> SLOPE = c(seq(0.2,0.6, .1),seq(0.2, -0.1, -0.1))) > >>>> > >>>> In this example, row 9 has a negative measurement for Chamber 5, > >>>> so I > >>>> would like to delete row 6, which is the same Chamber on the same > >>>> day, but not row 3, which is the same chamber on a different > >>>> day. In > >>>> the full dataframe, there are, of course, many more days. > >>>> > >>>> Is there a handy R way to do this? > >>>> > >>>> Thank you for the assistance. > >>>> > >>>> Toby > >>> > >>> > >>> > >>> Not fully tested, but here is one possibility: > >>> > >>>> toy > >>> CH DAY SLOPE > >>> 1 3 4 0.2 > >>> 2 4 4 0.3 > >>> 3 5 4 0.4 > >>> 4 3 4 0.5 > >>> 5 4 4 0.6 > >>> 6 5 5 0.2 > >>> 7 3 5 0.1 > >>> 8 4 5 0.0 > >>> 9 5 5 -0.1 > >>> > >>> > >>>> subset(toy, ave(SLOPE, CH, DAY, FUN = function(x) any(x < 0)) == 0) > >>> CH D
[R] syntax for batching rbind process
Dear helpeRs,
I am attempting to read in a series of csv files so I can bind them
into one large dataframe. I have written the following script:
test <- list.files(".", pattern = "csv") #lline 1
imp <- list()#line 2
for (i in 1:length(test)) {#line 3
imp[i] <- read.csv(test[i]) #line 4
}#line 5
works <- do.call(rbind, imp) #line 6
write.csv(works, "allmet.csv") #line 7
This script has an error at line 4. imp[i] contains only the value
of the first element of test[i]; in other words, every element of
imp[i] equals test[i] [1,1]. Otherwise, the script works. Could
someone please enlighten me as to the correct syntax for line 4?
Thank you in advance.
Toby
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] syntax for batching rbind process
Thank you for the suggestions for the more efficient code. The
problem remains, however, that the final dataframe does not contain
the correct values. So, in the case of the code you suggested,
imp <- lapply(test, read.csv)
do.call(rbind, imp)
imp does contain all the data from each dataframe, and the data from
each csv can be accessed with a single bracket index, but the do.call
does not work, possibly because rbind doesn't work on a list??
Any additional suggestions will be happily tested. I'm still
figuring out how to create a functionally equivalent toy example.
Thank you.
Toby
On 18 Aug 2010 at 13:49, Erik Iverson wrote:
>
>
> Toby Gass wrote:
> > Dear helpeRs,
> >
> > I am attempting to read in a series of csv files so I can bind them
> > into one large dataframe. I have written the following script:
> >
> > test <- list.files(".", pattern = "csv") #lline 1
> > imp <- list()#line 2
> > for (i in 1:length(test)) {#line 3
> > imp[i] <- read.csv(test[i]) #line 4
> > }#line 5
> > works <- do.call(rbind, imp) #line 6
> > write.csv(works, "allmet.csv") #line 7
> >
> > This script has an error at line 4. imp[i] contains only the value
> > of the first element of test[i]; in other words, every element of
> > imp[i] equals test[i] [1,1]. Otherwise, the script works. Could
> > someone please enlighten me as to the correct syntax for line 4?
>
> Just use lapply instead of a loop, easier syntax.
>
> imp <- lapply(test, read.csv)
> do.call(rbind, imp)
>
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] syntax for batching rbind process
It works perfectly now. Thank you all. Toby On 18 Aug 2010 at 15:04, Erik Iverson wrote: > > > Toby Gass wrote: > > Thank you for the suggestions for the more efficient code. The > > problem remains, however, that the final dataframe does not contain > > the correct values. So, in the case of the code you suggested, > > > > imp <- lapply(test, read.csv) > > do.call(rbind, imp) > > > > imp does contain all the data from each dataframe, and the data from > > each csv can be accessed with a single bracket index, but the do.call > > does not work, possibly because rbind doesn't work on a list?? > > > > I think that syntax looks OK. What error are you actually getting? > > > Any additional suggestions will be happily tested. I'm still > > figuring out how to create a functionally equivalent toy example. > > Toy example that can replicate error would be perfect. You can always > just ?dput a subset of imp, say the first few elements if there are > many. > __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] interpreting date-related error message
Hello, helpeRs, I have a vector of numbers from 1-365 (days of the year) that I would like to convert to a date. There are no NA's and no missing values. I did not insert leading zero's for numbers less than 100. Using the syntax: dat$doy.1 <- as.numeric(format(dat$doy, "%j" )) I get the following error message: Error in prettyNum(.Internal(format(x, trim, digits, nsmall, width, 3L, : invalid 'trim' argument .What is the error message telling me? (Windows OS and R 2.11.1) Thank you. Toby __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] interpreting date-related error message
Thanks to everyone for the explanations. Toby On 27 Aug 2010 at 12:46, Phil Spector wrote: > Toby - > Since dat$doy is just a number, the default S3 method > for format is used, where the second argument is the trim > parameter. I suspect you are confusing format (which is for > output) with strptime (which is for input). > For example, > > strptime(dat$doy,'%j') > > will assume that the dates are in the current year, and > return a POSIXlt object. Alternatively, you could pass > an origin to as.Date: > > as.Date(dat$doy,origin='2009-12-31') > > to get a similar Date object. > > - Phil > > > On Fri, 27 Aug 2010, Toby Gass wrote: > > > Hello, helpeRs, > > > > I have a vector of numbers from 1-365 (days of the year) that I would > > like to convert to a date. There are no NA's and no missing values. > > I did not insert leading zero's for numbers less than 100. > > > > Using the syntax: > > > > dat$doy.1 <- as.numeric(format(dat$doy, "%j" )) > > > > I get the following error message: > > > > Error in prettyNum(.Internal(format(x, trim, digits, nsmall, width, > > 3L, : > > invalid 'trim' argument > > > > .What is the error message telling me? > > (Windows OS and R 2.11.1) > > > > > > Thank you. > > > > Toby > > > > __ > > R-help@r-project.org mailing list > > https://stat.ethz.ch/mailman/listinfo/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] leap year and order function
Dear All,
I've always used this code:
year=c(1948:1953,2000,2100,2200,2300)
numdays=ifelse((year%%4==0 & year%%100!=0) | year%%400==0,366,365)
> numdays
[1] 366 365 365 365 366 365 366 365 365 365
Toby
From: r-help-boun...@r-project.org [r-help-boun...@r-project.org] On Behalf Of
Dennis Murphy [djmu...@gmail.com]
Sent: 31 January 2011 11:51
To: bill.venab...@csiro.au
Cc: R-help@r-project.org
Subject: Re: [R] leap year and order function
Hi:
There is an additional proviso:
If the year is divisible by 100, then the next test is whether it's
divisible by 400. If both conditions are met, then the century year is a
leap year. Hence, 2000 was a leap year because it is divisible by 400, but
2100, 2200 and 2300 will not be. For Perl and Ruby code, see
http://feb29.leapyearday.com/freecode.html
The 'correct' code is in a link associated with the Note below the first
paragraph.
HTH,
Dennis
On Sun, Jan 30, 2011 at 9:16 PM, wrote:
> > yearLength <- function(year) 365 + (year %% 4 == 0)
>
>
> > yearLength(1948:2010)
> [1] 366 365 365 365 366 365 365 365 366 365 365 365 366 365 365 365 366
> 365 365 365 366
> [22] 365 365 365 366 365 365 365 366 365 365 365 366 365 365 365 366 365
> 365 365 366 365
> [43] 365 365 366 365 365 365 366 365 365 365 366 365 365 365 366 365 365
> 365 366 365 365
> >
>
> -Original Message-
> From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org]
> On Behalf Of Bobby Lee
> Sent: Monday, 31 January 2011 2:23 PM
> To: R-help@r-project.org
> Subject: [R] leap year and order function
>
> im trying to write a for loop so that every leap year, the number of days
> becomes to 366 instead of 365. could someone help me out?
> and also, this set of data has 99.99s I set all the 99.99 ==NA.
> however, when im doing the order function to find the max value of that
> year, it still reads 99.99 as the max value.
> Thank you very much
>
> maxday <- matrix(NA,63,5)
>
> for (a in 1948:2010){
>
> maxday[,1]<-1948:2010
> yearly<-na.omit(dat.mat[dat.mat[,1]==a,])
>
> maxday[a-1947,2]<-yearly[order(yearly[,4])[*365*],2]
> maxday[a-1947,3]<-yearly[order(yearly[,4])[*365*],3]
>
>
> maxday[63,2]<-yearly[order(yearly[,4])[127],2]
> maxday[63,3]<-yearly[order(yearly[,4])[127],3]
> maxday[a-1947,4]<-max(yearly[,4])
> maxday[,5]<-len[,2]
> }
>
>[[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] how to learn more from R
Hi Jie Tang, You could try my "Friendly Beginners' R Course". Online resource at: http://www.R-project.org/ → Documentation/Other (left-hand panel) → Contributed Documentation (middle of screen) → (scroll down a little) which is a quick intro to what R can do in 14 pages with lots of links for further information. Toby Marthews From: r-help-boun...@r-project.org [r-help-boun...@r-project.org] On Behalf Of Jie TANG [totang...@gmail.com] Sent: 04 February 2011 06:42 To: r-help@r-project.org Subject: [R] how to learn more from R hi ,R USERS: I am a new R user. I have study some basic R syntax and get some interesting results. Now I want to know how much R can do .where can i get some examples or demo resource to study advanced function of R ? -- thank your TANG Jie Email: totang...@gmail.com Shanghai Typhoon Institute,China [[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] add error bars in a plot
Hi Maria, Instead of using errbar(), you can use the segments command which may be more flexible. Following the example on p.80 (the errbar command) of http://cran.ma.ic.ac.uk/web/packages/Hmisc/Hmisc.pdf set.seed(1) x=1:10 y=x+rnorm(10) delta=runif(10) barwidth=0.015 x11();plot(x,y) segments(x,y-delta,x,y+delta) segments(x-barwidth,y+delta,x+barwidth,y+delta) segments(x-barwidth,y-delta,x+barwidth,y-delta) Best, Toby From: r-help-boun...@r-project.org [r-help-boun...@r-project.org] On Behalf Of Gerrit Eichner [gerrit.eich...@math.uni-giessen.de] Sent: 09 February 2011 07:45 To: Lathouri, Maria Cc: r-help@r-project.org Subject: Re: [R] add error bars in a plot Hello, Maria, take a look at ?errbar in the package Hmisc. Hth -- Gerrit > Dear all > > I have a dataset of how metal concentrations change through time. I have > made a plot of date versus metal concentration. However I want to add > error bars in the plot. > > Could you help me? > > Thanks > Maria __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Linear regressions: producing multiple outputs
Hi RTSlider, I suspect you rather need to use the lme command (or perhaps glmmPQL or lmer) because you have a random predictor? lme(fixed=LeafLength~AirTemp*SnowFreeDate,random=~1|Species) See http://socserv.mcmaster.ca/jfox/Books/Companion-1E/appendix-mixed-models.pdf for a tutorial on lme. Toby From: r-help-boun...@r-project.org [r-help-boun...@r-project.org] On Behalf Of RTSlider [rob.t.sli...@gmail.com] Sent: 16 February 2011 18:13 To: r-help@r-project.org Subject: [R] Linear regressions: producing multiple outputs Hello all, I’m running simple linear regressions on multiple species of plants, comparing abiotic factor X against plant trait Y (e.g. Species1: leaf length vs air temperature). Ideally, what I’m looking for is an output giving me the R2, p value, coefficient, and Y intercept for each regression. Something like the example below: Species1: leaf length vs air temperature R2 = 0.10 p = 0.50m = 5.23b = 12.2 Species2: leaf length vs air temperature R2 = 0.10 p = 0.50m = 5.23b = 12.2 Species1: leaf length vs snow-free date R2 = 0.10 p = 0.50m = 5.23b = 12.2 Species2: leaf length vs snow-free date R2 = 0.10 p = 0.50m = 5.23b = 12.2 I currently have my data in this form: Species LeafLength AirTemp.SnowFreeDate Species11.1 20 160 Species24.5 20 160 Species35.4 20 160 And thought I could try this formula: lm(formula = LeafLength~AirTemp, SnowFreeDate | Species) But R is not a fan of it. Is there a way to do this (or get something remotely close to this)? I realize the output will probably be a bit messier than this, but what I’m really looking to avoid is running these regressions individually. -- View this message in context: http://r.789695.n4.nabble.com/Linear-regressions-producing-multiple-outputs-tp3309411p3309411.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] barplot with errorbars
If you google "barplot with error bars" you immediately find http://svitsrv25.epfl.ch/R-doc/library/prada/html/barploterrbar.html . Toby. From: r-help-boun...@r-project.org [r-help-boun...@r-project.org] On Behalf Of Lathouri, Maria [m.lathour...@imperial.ac.uk] Sent: 17 February 2011 16:00 To: r-help@r-project.org Subject: [R] barplot with errorbars Dear all I have six variables of the average metal concentrations Var1 4.77 Var2 23.5 Var3 5.2 Var4 12.3 Var5 42.1 Var6 121.2 I want to plot them as a barplot with error bars. Could you help me? Cheers Maria [[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] GSOC-R 2012
With John Nash, I am an administrator for the R project's contribution to the Google Summer of Code 2012. This is a program where Google provides $5000 to a student to write code for an open source project over a 3 month period in the summer. We would like to invite anyone interested in being either a mentor or a student intern to join our dedicated email list http://groups.google.com/group/gsoc-r and to check out the list of suggested projects http://rwiki.sciviews.org/doku.php?id=developers:projects:gsoc2012 Please add to the table of projects on the RWiki if you have any ideas. Thanks very much! [[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] glmulti fails because of rJava
Dear R,
The glmulti package no longer loads through the library() command, apparently
because of a problem with rJava.
I have today reinstalled R from scratch (updated to v2.14.0) and reinstalled
all packages from scratch and updated them all too. The problem is the same as
I found on v2.13.2. See session below for the error. I tried
install.packages(rJava) as advised by the error report but it didn't help (see
session below).
Any advice would be much appreciated: I can now not use glmulti at all and I
can't see how to solve this!
Thank you very much R Help!
Best,
Toby
R version 2.14.0 (2011-10-31)
Copyright (C) 2011 The R Foundation for Statistical Computing
ISBN 3-900051-07-0
Platform: x86_64-pc-mingw32/x64 (64-bit)
R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.
Natural language support but running in an English locale
R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.
Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.
[Previously saved workspace restored]
> library("glmulti")
Loading required package: rJava
Error : .onLoad failed in loadNamespace() for 'rJava', details:
call: stop("No CurrentVersion entry in '", key, "'! Try re-installing Java
and make sure R and Java have matching architectures.")
error: object 'key' not found
Error: package ‘rJava’ could not be loaded
> install.packages("glmulti")
--- Please select a CRAN mirror for use in this session ---
trying URL
'http://www.stats.bris.ac.uk/R/bin/windows/contrib/2.14/glmulti_1.0.2.zip'
Content type 'application/zip' length 104179 bytes (101 Kb)
opened URL
downloaded 101 Kb
package ‘glmulti’ successfully unpacked and MD5 sums checked
The downloaded packages are in
C:\Users\TobyM\AppData\Local\Temp\RtmpUD3kYK\downloaded_packages
> library("glmulti")
Loading required package: rJava
Error : .onLoad failed in loadNamespace() for 'rJava', details:
call: stop("No CurrentVersion entry in '", key, "'! Try re-installing Java
and make sure R and Java have matching architectures.")
error: object 'key' not found
Error: package ‘rJava’ could not be loaded
> install.packages("rJava")
trying URL
'http://www.stats.bris.ac.uk/R/bin/windows/contrib/2.14/rJava_0.9-2.zip'
Content type 'application/zip' length 676445 bytes (660 Kb)
opened URL
downloaded 660 Kb
package ‘rJava’ successfully unpacked and MD5 sums checked
The downloaded packages are in
C:\Users\TobyM\AppData\Local\Temp\RtmpUD3kYK\downloaded_packages
> library("glmulti")
Loading required package: rJava
Error : .onLoad failed in loadNamespace() for 'rJava', details:
call: stop("No CurrentVersion entry in '", key, "'! Try re-installing Java
and make sure R and Java have matching architectures.")
error: object 'key' not found
Error: package ‘rJava’ could not be loaded
>
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
[R] Constructing groupedData objects in nlme - a little problem
Dear R-help,
I am trying to create groupedData objects using the nlme library. I'm
missing something basic, I know:
Here is the first example in ch.1 of Pinheiro & Bates (2000):
library(nlme)
x2=Rail$travel;x1=Rail$Rail;eg1=data.frame(x1,x2);eg1gd=Rail
print(eg1gd)
x11();print(plot(eg1gd))
femodel=lm(x2~x1-1,data=eg1gd)
print(femodel$coefficients)
Result:
x12 x15 x11 x16 x13 x14
31.7 50.0 54.0 82.7 84.7 96.0
...which works fine. This uses a built-in groupedData object called "Rail"
that is part of the nlme library.
I am trying to 'recreate' this groupedData object. Here's what I've done:
x1=c(1,1,1,2,2,2,3,3,3,4,4,4,5,5,5,6,6,6);x2=c(55,53,54,26,37,32,78,91,85,92,100,96,49,51,50,80,85,83)
eg1=data.frame(x1,x2);colnames(eg1)=c("Rail","travel");eg1gd=groupedData(travel~1|Rail,data=eg1)
print(eg1gd)
x11();print(plot(eg1gd))
femodel=lm(x2~x1-1,data=eg1gd)
print(femodel$coefficients)
Result:
x1
16.49817
...but, as you can see, the coefficients I get at the end this time are
completely different and I don't know why. Somehow, I am not creating the
structure properly even though the formula and data values are all
correct.
Can anyone help? I've looked at the ?groupedData man page, but it has no
solution to this.
Thanks very much for any advice,
Toby
Pinheiro JC & Bates DM (2000). Mixed-Effects Models in S and S-PLUS (1st
ed.). Springer, New York.
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Word wrapping for character objects (WINDOWS R ONLY)
Can anybody help me with this problem? ** ONLY WINDOWS R - PROBLEM DOESN'T
OCCUR ON LINUX **
I want to print a long character to screen:
> getOption("width")
[1] 60
> z=(1:20)/10#z is a vector of length between 20 and 30 (depending on
user options) containing lengths in mm (i.e. each element is 1-5
characters long)
> str1=paste("The depths chosen are (",toString(z),") mm, and more text
...\n")
> cat(str1)
The depths chosen are ( 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1,
1.1, 1.2, 1.3, 1.$
>
The problem is that on R for Windows the string is cropped by the window
size (hence the "$"). On R for Linux, this doesn't happen and the text is
"word wrapped" (the default for the shell, I guess):
> cat(str1)
The depths chosen are ( 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1,
1.1, 1.2, 1.3,
1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2 ) mm, and more text ...
>
I can't find any option for "word wrapping" in the cat command (fill=TRUE
has no effect). I also checked the menu Edit -> GUI preferences..., but
there doesn't seem to be a "Word Wrap" option there either.
How do I get word wrapping like this in Windows? Perhaps the attached
screenshots clarify this question.
THANKS FOR ANY HELP!
Toby Marthews
Previous relevant posts:
- The post from 2006 about Screen Wrapping
(http://tolstoy.newcastle.edu.au/R/help/06/05/26673.html) which Brian
Ripley answered was about controlling how long vectors are cropped to the
screen. Unfortunately, the width option in options() does not affect
character objects, so I can't use that control here.
- I sent the same question to [EMAIL PROTECTED] in Oct 2007, but
noone there could help me with it.
- Try the following command on Windows R with a small window
(getOption("width")<117) and a large window (getOption("width")>117) and
you'll see you get extra nonexistent options in the menu:
a=c("Jan","Feb","Mar","Apr","May","Jun","Jul","Aug","Sep","Oct","Nov","Dec");menu(a)
I guess this is a related problem?__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] controlling location of labels in axis()
Hi Andrew,
Perhaps this example would help. You can add in spaces to the mtext text
to move the text sideways.
par(mai=c(0.5,0.5,0.5,0.5),oma=c(2,2,2,2)) #mai units are INCHES, oma
units are LINES
plot(runif(50),xlab="xlab",ylab="ylab",bty="l") #n.b. these labels don't
appear
mtext("First inner x axis label",side=1)
mtext("First inner y axis label",side=2)
mtext("Second inner x axis label",side=3)
mtext("Second inner y axis label",side=4)
mtext("First outer x axis label",side=1,outer=TRUE)
mtext("First outer y axis label",side=2,outer=TRUE)
mtext("Second outer x axis label",side=3,outer=TRUE)
mtext("Second outer y axis label",side=4,outer=TRUE)
Cheers,
Toby Marthews
Le Jeu 12 juin 2008 15:32, Andrew Yee a écrit :
> Here's a naive question about axis()
>
> How do you control the location of the labels with the axis() command?
>
> In the following example:
>
> foo <- data.frame(plot.x=seq(1:3), plot.y=seq(4:6))
> plot(foo$plot.x, foo$plot.y, type='n', axes=FALSE)
> points(foo$plot.x, foo$plot.y)
> axis(1, at=foo$plot.x, labels=foo$plot.x)
>
> I'd like to be able the control the y location of the labels (i.e. move
up
> or down in the graph).
>
> Thanks,
> Andrew
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Level Plot and Scale of Colorkey
Try colscaledivs=100#colscaledivs=15 here is the R default levelplot(z ~ x * y, g,xlab="x co-ordinate",ylab="y co-ordinate",colorkey=TRUE,at=seq(from=-0.01,to=0.25,length=colscaledivs),col.regions=(col=gray((0:colscaledivs)/colscaledivs))) Toby Marthews Le Ven 13 juin 2008 18:50, emma hartnett a écrit : > I am drawing level plots but I would like to specify the range of the > colorkey, I am not having any success figuring this out so any help would > be greatly appreciated! > > Here is an example of what I am trying to do: > > disp<-1 > > x <- seq(1, 10,by=1) > y <- seq(1,10,by=1) > g <- expand.grid(x = x, y = y) > g$z <- 1/exp((abs(g$x-5)+abs(g$y-5))*disp) > g$z<-g$z/sum(g$z) > > levelplot(z ~ x * y, g,xlab="x co-ordinate", ylab="y co-ordinate" > ,colorkey=TRUE,col.regions=(col=gray((0:32)/32))) > > I would like to enforce the number of divisions on the colorkey scale and > the size – so for example from 0 to 0.1 in increments of 0.02 (just as an > example). > > I apologize if this is an obvious question but I have read the > documentation and scoured the archives and cannot figure it out. > > > > > __ > can.html > > __ > R-help@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide > http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] round(1.5) = round(2.5) = 2?
Hi Markus, The R function round() uses the round-to-even method, which is explained on http://en.wikipedia.org/wiki/Rounding#Round-to-even_method. If you would like instead "traditional rounding" then you should add 0.5 and take the integer part, as is suggested in the examples on ?round, e.g. x=c(1.49,1.50,1.51,2.49,2.50,2.51,3.49,3.50,3.51);r2e=round(x);trad=floor(x+0.5);data.frame(x,r2e,trad) The reason for the IEEE standard is to do with signal processing and the bias introduced by traditional rounding if you have a lot of data points whose decimal expansion ends in 5. Personally, I am a mathematician and satistician and I find that in 99% of cases traditional rounding is what is required (basically always except some very specific examples involving very large sets of data) so for R to have put the IEEE standard as the default (rather than, say, an option) is a bit odd. However, R's benefits and advantages by far outweigh its little oddities, as I presume you know since you are using it. Effectively, I never use the round() command and always calculate using the floor function. Toby Marthews Le Dim 15 juin 2008 11:26, Markus Didion a écrit : > Dear R-users > > with a bit of grief I had to repeat an extensive analysis because I > did not suspect (and therefore did not read the documentation) that > round was implemented as "for rounding off a 5, the IEC 60559 standard > is expected to be used, 'go to the even digit'", resulting in > round(1.5) = 2 > round (2.5) = 2. > > As a non-mathematician I am both puzzled and intrigued by this rule as > it is against what I have learned in my math courses, i.e. > round(1.5) = 2 > round (2.5) = 3. > > I would like to understand the reason behind this rule. > > Thanks for your comments. > > Markus > > -- > > Markus Didion > > Wald?kologie Forest Ecology > Inst. f. Terrestrische OekosystemeInst. of Terrestrial Ecosystems > Departement Umweltwissenschaften Dept. of Environmental Sciences > Eidg. Technische Hochschule Swiss Fed. Inst. of Technology > ETH-Zentrum CHN G78 ETH-Zentrum CHN G78 > Universit?tstr. 22 > Universitaetstr. 22 > CH-8092 Z?richCH-8092 Zurich > Schweiz > Switzerland > > Tel +41 (0)44 632 5629Fax +41 (0)44 632 1358 > Email [EMAIL PROTECTED] > homepage: http://www.fe.ethz.ch/people/didionm __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] How can I stop strwrap removing escape characters? (and solution to the screen wrapping question from 11/6)
After the helpful reply from Greg Snow (below), I've written a function
printtoscreen which does the screen wrapping I need (although not in a
very elegant way). This works fine for strings without any escape
characters in them, e.g.
> printtoscreen=function(str) {
+ str2=strwrap(str,width=getOption("width"))
+ if (length(str2)>1)
{str2[2:length(str2)]=paste("\b\b",str2[2:length(str2)],sep="")} #to
remove the ", "s toString puts in
+ str2=paste(str2,"\n",sep="")
+ cat(str2)
+ }
> getOption("width")
[1] 59
> cat("Measured lengths are (",toString(1:20/100),") mm.\n")
Measured lengths are ( 0.01, 0.02, 0.03, 0.04, 0.05, 0.06,$
> printtoscreen(paste("Measured lengths are (",toString(1:20/100),")
mm.\n",sep=""))
Measured lengths are (0.01, 0.02, 0.03, 0.04, 0.05, 0.06,
0.07, 0.08, 0.09, 0.1, 0.11, 0.12, 0.13, 0.14, 0.15, 0.16,
0.17, 0.18, 0.19, 0.2) mm.
>
HOWEVER, if you use any escape characters in your string, it fails because
strwrap removes all of them:
>
> cat("\n\tHello.\n")
Hello.
> printtoscreen("\n\tHello.\n")
Hello.
>
Is there any way to prevent strwrap doing this? I haven't found any
options for this.
Thanks very much (and thanks to Greg Snow for the tip about strwrap before).
Toby Marthews
==
Le Mer 11 juin 2008 18:00, Greg Snow a écrit :
> You could try passing your character string to the strwrap function
first,
> then use cat on the result.
>
> --
> Gregory (Greg) L. Snow Ph.D.
> Statistical Data Center
> Intermountain Healthcare
> [EMAIL PROTECTED]
> (801) 408-8111
>
>
>
>> -Original Message-
>> From: [EMAIL PROTECTED]
>> [mailto:[EMAIL PROTECTED] On Behalf Of Toby Marthews Sent:
Wednesday, June 11, 2008 6:52 AM
>> To: r-help@r-project.org
>> Subject: [R] Word wrapping for character objects (WINDOWS R ONLY) Can
anybody help me with this problem? ** ONLY WINDOWS R - PROBLEM
DOESN'T OCCUR ON LINUX **
>> I want to print a long character to screen:
>> > getOption("width")
>> [1] 60
>> > z=(1:20)/10#z is a vector of length between 20 and 30 (depending
on user options) containing lengths in mm (i.e. each element is 1-5
characters long)
>> > str1=paste("The depths chosen are (",toString(z),") mm, and more text
...\n")
>> > cat(str1)
>> The depths chosen are ( 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1,
1.1, 1.2, 1.3, 1.$
>> >
>> The problem is that on R for Windows the string is cropped by the
window size (hence the "$"). On R for Linux, this doesn't happen and the
text is "word wrapped" (the default for the shell, I guess):
>> > cat(str1)
>> The depths chosen are ( 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1,
1.1, 1.2, 1.3,
>> 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2 ) mm, and more text ...
>> >
>> I can't find any option for "word wrapping" in the cat command
(fill=TRUE has no effect). I also checked the menu Edit -> GUI
preferences..., but there doesn't seem to be a "Word Wrap" option there
either.
>> How do I get word wrapping like this in Windows?
>> THANKS FOR ANY HELP!
>> Toby Marthews
>> Previous relevant posts:
>> - The post from 2006 about Screen Wrapping
(http://tolstoy.newcastle.edu.au/R/help/06/05/26673.html) which Brian
Ripley answered was about controlling how long vectors are cropped to the
screen. Unfortunately, the width option in options() does not
affect character objects, so I can't use that control here.
>> - I sent the same question to [EMAIL PROTECTED] in Oct 2007,
but noone there could help me with it.
>> - Try the following command on Windows R with a small window
(getOption("width")<117) and a large window (getOption("width")>117) and
you'll see you get extra nonexistent options in the menu:
>> a=c("Jan","Feb","Mar","Apr","May","Jun","Jul","Aug","Sep","Oct","Nov","Dec");menu(a)
I guess this is a related problem?
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Calling functions
Hi,
Can't you just do source("Xtabs.R")? That will load in the definition.
Alternatively, instead of saving the R program, save the workspace on your
network (e.g. Xtabs.RData), which will contain the function definition,
and arrange a link so that R always starts with that workspace loaded?
Toby Marthews
Le Mar 17 juin 2008 13:32, Michael Pearmain a écrit :
> Another newbie question.
>
> I've written a function and saved the file as Xtabs.R, in a central place
> on
> a network so others will be able ot use the function,
> My question is how do i call this function?
>
> I've tried to chance the working directory, and tried to load it via;
>> library(Xtabs, lib.loc="//filer/common/technical/surveys/R_test")
>
> but neither seem to work? the function inside is called CrossTable.
>
> Mant thanks in advance
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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
Hi there, Simplifying your example a bit, we have: > a=data.frame(r=c(0,0,1,1),g=c(0,1,0,1)) > transform(a,R=(r>0 && g>0)) r g R 1 0 0 FALSE 2 0 1 FALSE 3 1 0 FALSE 4 1 1 FALSE > transform(a,R=(r>0 & g>0)) r g R 1 0 0 FALSE 2 0 1 FALSE 3 1 0 FALSE 4 1 1 TRUE > ...but actually, in the && case the R column values are not calculated from the corresponding values of r and g: the FALSE is effectively calculated in row 1 and then copied down. Compare with this: > a=data.frame(r=c(1,1,0,0),g=c(1,0,1,0)) > transform(a,R=(r>0 && g>0)) r gR 1 1 1 TRUE 2 1 0 TRUE 3 0 1 TRUE 4 0 0 TRUE > transform(a,R=(r>0 & g>0)) r g R 1 1 1 TRUE 2 1 0 FALSE 3 0 1 FALSE 4 0 0 FALSE > This is actually explained in the ?!base man page: "& and && indicate logical AND and | and || indicate logical OR. The shorter form performs elementwise comparisons in much the same way as arithmetic operators. The longer form evaluates left to right examining only the first element of each vector. Evaluation proceeds only until the result is determined." Hope this helps, Toby Marthews Le Mer 18 juin 2008 15:10, Christos Argyropoulos a écrit : > > Hi, > > I noticed whether some one could explain why "&" and "&&" behave > differently in data frame transformations. > > Consider the following : > > a<-data.frame(r=c(0,0,2,3),g=c(0,2,0,2.1)) > > Then: > >> transform(a,R=ifelse(r>0 && g> 0,log(r/g),NA)) > > r g R > 1 0 0.0 NA > 2 0 2.0 NA > 3 2 0.0 NA > 4 3 2.1 NA > > but > >> transform(a,R=ifelse(r>0 & g> 0,log(r/g),NA)) > r g R > 1 0 0.0NA > 2 0 2.0NA > 3 2 0.0NA > 4 3 2.1 0.3566749 > > > If my understanding of the differences between "&" and "&&" and how > 'transform' works are accurate, both statements should produce the same > output. > > > I got the same behaviour in Windows XP Pro 32-bit (running R v 2.7) and > Ubuntu Hardy (running the same version of R). > > > Thanks > > Christos Argyropoulos > > University of Pittsburgh Medical Center __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 message: Bad value
Dear R help, In the middle of my session, I started receiving the error message "bad value" regardless of what I entered. I had to close R and restart it. When the error occurred, I was creating basic lattice plots, with no other add-on packages running. I am running R through John Fox's XEmacs interface, XEmacs version 21.4.20 on a Windows XP Pro machine with 2 GB of RAM. I have had no other problems. I've found 2 similar reports in the R-help archive, neither of which was resolved; one of the earlier posters was running Linux. platform i386-pc-mingw32 arch i386 os mingw32 system i386, mingw32 status major 2 minor 6.2 year 2008 month 02 day 08 svn rev 44383 language R version.string R version 2.6.2 (2008-02-08) Previous reports are found here: http://tolstoy.newcastle.edu.au/R/e4/help/08/01/0686.html http://tolstoy.newcastle.edu.au/R/e2/help/07/06/19095.html I do not need a response, unless there is a preventative measure, but thought I would bring this to your attention. Thank you, Toby Toby Gass Graduate Degree Program in Ecology Department of Forest, Rangeland, and Watershed Stewardship Warner College of Natural Resources Colorado State University Fort Collins, CO 80523 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 message: Bad value
Thank you for the response. I was doing a basic operation with a small dataset and am (fortunately for me) unable to reproduce the error. After restarting R, without shutting down anything else, there was no error. I thought I would post, however, since 2 others had reported the same problem in the past year. Toby - Original Message - From: "Duncan Murdoch" <[EMAIL PROTECTED]> To: "Toby Gass" <[EMAIL PROTECTED]> Cc: Sent: Tuesday, June 24, 2008 6:27 PM Subject: Re: [R] Error message: Bad value | On 24/06/2008 3:12 PM, Toby Gass wrote: | > Dear R help, | > | > In the middle of my session, I started receiving the error | > message "bad value" | > regardless of what I entered. I had to close R and restart | > it. When the error | > occurred, I was creating basic lattice plots, with no other | > add-on packages running. | > I am running R through John Fox's XEmacs interface, XEmacs | > version | > 21.4.20 on a Windows XP Pro machine with 2 GB of RAM. I | > have had | > no other problems. I've found 2 similar reports in the | > R-help | > archive, | > neither of which was resolved; one of the earlier posters | > was running Linux. | | That message is probably being printed in response to memory corruption. | | > platform i386-pc-mingw32 | > arch i386 | > os mingw32 | > system i386, mingw32 | > status | > major 2 | > minor 6.2 | | That version is a little old; 2.7.1 was just released. If you can make | the error reproducible in a current version (or even in 2.6.2), it would | be very helpful to post the recipe, and it will probably get fixed. | | Duncan Murdoch | | | > year 2008 | > month 02 | > day 08 | > svn rev 44383 | > language R | > version.string R version 2.6.2 (2008-02-08) | > | > Previous reports are found here: | > http://tolstoy.newcastle.edu.au/R/e4/help/08/01/0686.html | > http://tolstoy.newcastle.edu.au/R/e2/help/07/06/19095.html | > | > I do not need a response, unless there is a preventative | > measure, | > but thought I would bring this to your attention. | > | > Thank you, | > | > Toby | > | > Toby Gass | > Graduate Degree Program in Ecology | > Department of Forest, Rangeland, and Watershed Stewardship | > Warner College of Natural Resources | > Colorado State University | > Fort Collins, CO 80523 | > | > __ | > R-help@r-project.org mailing list | > https://stat.ethz.ch/mailman/listinfo/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] nested time series data with measurement error
Dear R helpers,
I am trying to fit a model with the main objective of
assessing differences, rather than predicting.
The treatment was applied to half of the subjects after 9
months of measurement.
Then 9 more months of data were collected. The variable
"month" has measurement error
due to the complexities of data collection.
The dependent variable is numeric and continuous.
This would be a simple difference between means were it not
for the repeated
measurements and measurement error.
I am wondering what the R gurus suggest as to which R
function might best model
these data. I have been looking at tsls (two-stage least
squares) in the sem package
and pls (partial least squares). I am not sure about their
compatibility with repeated
measurements and time series. The response curve is also
likely to be non-linear.
The following script approximates the data, including sample
size. In the real data, there are 18 months spread out over
a
36 month period.
#generate data
tree<- gl(6,18,label = paste("tree",1:6))
month <- gl(18,1,length = 108,label = paste("month",1:18),
ordered = TRUE)
trtmt <- gl(2, 54, length = 108, label = paste("trt",1:2))
pre.post <- gl(2,9,length = 108, label = c("pre","post"))
response <- runif(108, min = -28, max = -25)
help <- data.frame(tree,month,trtmt,pre.post, response)
Thank you in advance for your assistance.
Toby Gass
Graduate Degree Program in Ecology
Department of Forest, Rangeland, and Watershed Stewardship
Warner College of Natural Resources
Colorado State University
Fort Collins, CO 80523
email: [EMAIL PROTECTED]
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
