[R] How to set up a Negative Exponential regression mixed effect model 'nlme' ?

2016-11-17 Thread Ronny Steen
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.

2015-11-22 Thread Ronny Steen
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

2013-08-12 Thread ronny
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

2013-01-20 Thread ronny
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

2012-09-13 Thread ronny
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

2011-11-15 Thread ronny
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

2011-07-19 Thread ronny
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

2011-07-19 Thread ronny
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'

2008-11-21 Thread Ronny Wölbing




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'

2008-11-20 Thread Ronny Wölbing

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?

2008-01-29 Thread Ronny Steen

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.