[R] How to set up a Negative Exponential regression mixed effect model 'nlme' ?
Hi, I need help with fitting a non-linear mixed effects model (nested random effect). I look at the relationship between 'wingbeat frequency' (beats per sec) and 'wing length' in hummingbirds, and I will find a model with the best fit. I have posted the question at R-SIG-mixed-models list, but no answer so far, therefore I give it a new try here. I manage to fit Negative Natural Exponential regression and Assymptotic (SSasymp-function) nonlinear regression, but I don't figure out how to fit Negative Exponential regression (similar to nls {stats} function 'negexp.SSival'). *My question:* Could I get some help adapting the Negative Exponential regression model (*'nls'*, without random effect) to a model with nested random effect (*'nlme' *)? Script and access to data: library(nlme) library(MASS) WBF <-read.csv(url("https://www.dropbox.com/s/hin8o27i1kmdloe/ Species2016.csv.csv?raw=1")) plot(WBF$WL,WBF$Beat_freq,type="p",ylab="WBF")# Make a plot #Negative Natural Exponential regression with nested random effect fm1<-lme(Beat_freq ~ exp(-WL), data = WBF, random = ~ 1|Species/ID, method = "ML") #Assymptotic (SSasymp-function) nonlinear regression with nested random effect fm2 <- nlme(Beat_freq ~ SSasymp(WL, Asym, R0, lrc), data = WBF, fixed = Asym + R0 + lrc ~ 1, random = Asym ~ 1|Species/ID, start = c(Asym = 27, R0 = 6000, lrc = 0.19)) #Negative Exponential regression - HELP NEEDED TO INCLUDE NESTED RANDOM EFFECT ('nlme') negexp<-selfStart(model = ~b0 + b1*exp(-x/th), initial = negexp.SSival,parameters=c("b0","b1","th"), template=function(x,b0,b1,th){}) fm3<-nls(Beat_freq~negexp(WL,B0,B1,theta), data=WBF,control=list(maxiter =5000),trace=T) Regards, Ron [[alternative HTML version deleted]] __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see https://stat.ethz.ch/mailman/listinfo/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] Converting time zones in R using metadata file information of video files, help needed.
Hi, I have video files (FAT) that are taken in a different timezone than my current location. The modification date/time in the metafila data of the video file shows the time video was taken, although in the current timezone of my computer, if I understand right. I wish to convert the date/time to the origin. The video was taken in London, Ontario Canada at 2015-06-21 07:53:28, when looking at the metadata of the file on my computer (timezone "Europe/Berlin") it says modification date "2015-06-22 01:53:28". Hence, there is a 6 hour difference between the two time-zones I use the script provided here: http://blog.revolutionanalytics.com/2009/06/converting-time-zones.html pb.txt <- "2015-06-22 01:53:28" #modification date as shown on my computer (timezone "Europe/Berlin") pb.date <- as.POSIXct(pb.txt, tz="Europe/London") pb.date <- as.POSIXct(pb.txt, tz="America/Toronto")#timezone of origin video camera location format(pb.date, tz=local.time,usetz=TRUE) #formats the mtime to data and time when the video was taken [1] "2015-06-22 07:53:28 CEST" # the time is correct but the date is wrong as it has added 6 hours rather than subtracting. I find working with different time-zones to be difficult, I hope I managed to formulate an understandable question. Regards, Kes [[alternative HTML version deleted]] __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see https://stat.ethz.ch/mailman/listinfo/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] degrees of freedom in custom contrasts ananlysis
Hi, I would like to understand why the residual standard error, and the degrees of freedom are changing when I define custom contrasts, which are not orthogonal. For example: y <- rnorm(40) x <- factor(rep(1:10,4)) summary(lm(y~x)) #standard model: Residual standard error: 1.103 on 30 degrees of freedom summary(lm(y~C(x,contr.sum))) #orthogonal contrasts: Residual standard error: 1.103 on 30 degrees of freedom summary(lm(y~C(x,c(rep(0.2,5),1,rep(0,4)),1))) #custom contrast: Residual standard error: 1.049 on 38 degrees of freedom In addition, is there a way to run the custom contrasts analysis without equal variance assumption? I expect the degrees of freedom to get lower in this case, as in ANOVA welch. Thanks! -- View this message in context: http://r.789695.n4.nabble.com/degrees-of-freedom-in-custom-contrasts-ananlysis-tp4673552.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] Lattice levelplot- remove unused levels per panel
Hi, I am using levelplot, and would like remove from each panel (condition) its unused x levels. e.g. Remove from panel vs=1 the cyl level=8. data(mtcars) levelplot(mpg~factor(cyl)*factor(gear)|factor(vs)) Thanks for your help, Ronny -- View this message in context: http://r.789695.n4.nabble.com/Lattice-levelplot-remove-unused-levels-per-panel-tp4656087.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] problematic p values for Dixon Q test
I am using R dixon.test in order to perform Q test for detection of outliers. The p-values for high Q values are consistent with the table in Rorabacher, D.B. (1991) http://www.flworkshop.com/sscs/dixon_test.pdf, but it seems that small Q-values get problematic p-values. For example: dixon.test(c(0.324,0.5,1,1.324),two.sided=T) Q = 0.324, p-value = 1 dixon.test(c(0.324,0.5,1,1.324),two.sided=F) Q = 0.324, p-value = 0.5 The p value I expect by simulation is ~0.83: norm=matrix(nrow=1000,ncol=4) for (i in 1:1000){ norm[i,]<-rnorm(4) q[i]<-dixon.test(norm[i,])$statistic p[i]<-dixon.test(norm[i,])$p.value>} hist(q) expected.p.val<- length(which(q>0.324))/1000 In particular I am concerened from the misleading result of dixon.test(c(1,1,2,2),two.sided=T) Q = 0, p-value < 2.2e-16 I would appreciate any help with this issue, -- View this message in context: http://r.789695.n4.nabble.com/problematic-p-values-for-Dixon-Q-test-tp4642983.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] Bootstrap values for hierarchical tree based on distaance matrix
I would like to get an hierarchical clustering tree with bootstrap values indicated on the nodes, as in pvclust. The problem is that I have only distance matrix instead of the raw data, required for pvclust. Is there a way to get it? fit1 <- hclust(dist) # an object of class '"dist" plot(fit1) # dendogram without p values library(pvclust) fit2 <- pvclust(raw.data, method.hclust="ward", method.dist="euclidean") plot(fit2) # dendogram with p values -- View this message in context: http://r.789695.n4.nabble.com/Bootstrap-values-for-hierarchical-tree-based-on-distaance-matrix-tp4042939p4042939.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Centering data frame by factor
Perfect! Made my day! -- View this message in context: http://r.789695.n4.nabble.com/Centering-data-frame-by-factor-tp3677609p3677665.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] Centering data frame by factor
Hi, I would like to center P1 and P2 of the following data frame by the factor "Experiment", i.e. substruct from each value the average of its experiment, and keep the original data structure, i.e. the experiment and the group of each value. RAW= data.frame("Experiment"=c(2,2,2,1,1,1),"Group"=c("A","A","B","A","A","B"),"P1"=c(10,12,14,5,3,4),"P2"=c(8,12,16,2,3,4)) Desired result: NORMALIZED= data.frame("Experiment"=c(2,2,2,1,1,1),"Group"=c("B","A","B","B","A","B"),"P1"=c(-2,0,2,1,-1,0),"P2"=c(-4,0,4,-1,0,1)) I tried using "by", but then I lose the original order, and the "Group" varaible. Can you help? > RAW Experiment Group P1 P2 2 A 10 8 2 A 12 12 2 B 14 16 1 A 5 2 1 A 3 3 1 B 4 4 NOT.OK<- within (RAW, {P1<-do.call(rbind,by(RAW$P1,RAW$Experiment,scale,scale=F))}) > NOT.OK Experiment Group P1 P2 2 A 1 8 2 A -1 12 2 B 0 16 1 A -2 2 1 A 0 3 1 B 2 4 -- View this message in context: http://r.789695.n4.nabble.com/Centering-data-frame-by-factor-tp3677609p3677609.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] How to convert S4 class slots into data.frame or how to assign variables of type 'Date'
Hi, I created a class (S4) with some slots like value, date, description (it's actually a financial transaction class). Now I need a method to convert this class forth and back into a single row data.frame, where every slots represents a column. This method looks at the moment like this: > setMethod("as.data.frame", "Transaction", function(x, row.names = NULL, optional = FALSE, ...){ slotnames <- slotNames(x) slotlist <- data.frame(rbind(1:length(slotnames))) names(slotlist) <- slotnames for(i in slotnames) { slotlist[1, i] <- slot(x, i) } return(slotlist) } ) This method doesn't require predetermined slotnames or types, which is important to me. The method works quite good but the problem is that I have slots of type 'Date' and this method doesn't preserve the type but converts it to numeric. A couple of tests showed that this is actually a problem of assigning values to data.frame column. Something like this: > slotlist$DayOfTransaction <- slot(Transaction, DateOfTransaction) would preserve the type of DateOfTransaction as 'Date'. But I don't see a way to use this assigning scheme in my method without using the actual slotnames and giving up a lot of flexibility. Do you have any suggestions? Is there maybe even a simple way to convert S4 slots into data.frames? Thanks in advance Ronny __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Does nobody has an idea? __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] How to convert S4 class slots into data.frame or how to assign variables of type 'Date'
Hi, I created a class (S4) with some slots like value, date, description (it's actually a financial transaction class). Now I need a method to convert this class forth and back into a single row data.frame, where every slots represents a column. This method looks at the moment like this: > setMethod("as.data.frame", "Transaction", function(x, row.names = NULL, optional = FALSE, ...){ slotnames <- slotNames(x) slotlist <- data.frame(rbind(1:length(slotnames))) names(slotlist) <- slotnames for(i in slotnames) { slotlist[1, i] <- slot(x, i) } return(slotlist) } ) This method doesn't require predetermined slotnames or types, which is important to me. The method works quite good but the problem is that I have slots of type 'Date' and this method doesn't preserve the type but converts it to numeric. A couple of tests showed that this is actually a problem of assigning values to data.frame column. Something like this: > slotlist$DayOfTransaction <- slot(Transaction, DateOfTransaction) would preserve the type of DateOfTransaction as 'Date'. But I don't see a way to use this assigning scheme in my method without using the actual slotnames and giving up a lot of flexibility. Do you have any suggestions? Is there maybe even a simple way to convert S4 slots into data.frames? Thanks in advance Ronny __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Guidance for reporting results from lme test?
Hello, I have sucsessfully used a linear mixed effect model, lme, (REML). The results are satisfactory, but I have problems with sorting out how to report the result in a scientific paper. Is there a genearal guidance for reporting lme results available (web page, book or article)? Best regards Ronny _ [[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.