On 7/20/2006 7:16 PM, Thomas Lumley wrote: > On Thu, 20 Jul 2006, [EMAIL PROTECTED] wrote: > >> Duncan, i could not find one seed which caused both corGaus and corExp >> to crash in R-patched. but, i found a seed that caused each to fail >> individually. Thanks for your help. > > > Running these under Valgrind they both show the same problem > > ==10878== Invalid read of size 4 > ==10878== at 0x1CCDB9EB: mult_mat (matrix.c:84) > ==10878== by 0x1CCD81B9: corStruct_recalc (corStruct.c:72) > ==10878== by 0x1CCDC06E: nlme_wtCorrAdj (nlme.c:152) > ==10878== by 0x1CCDC4F7: fit_nlme (nlme.c:347) > ==10878== by 0x80A518E: do_dotCode (dotcode.c:1777) > ==10878== by 0x80C45A9: Rf_eval (eval.c:444) > ==10878== by 0x80C600A: do_set (eval.c:1340) > ==10878== by 0x80C44A2: Rf_eval (eval.c:427) > ==10878== by 0x80C6097: do_begin (eval.c:1104) > ==10878== by 0x80C44A2: Rf_eval (eval.c:427) > ==10878== by 0x80C61AD: do_repeat (eval.c:1066) > ==10878== by 0x80C44A2: Rf_eval (eval.c:427) > ==10878== Address 0x1D4A60DC is 4 bytes after a block of size 6048 alloc'd > ==10878== at 0x1B909B71: calloc (vg_replace_malloc.c:175) > ==10878== by 0x80E09D5: R_chk_calloc (memory.c:2315) > ==10878== by 0x1CCDC415: fit_nlme (nlme.c:113) > ==10878== by 0x80A518E: do_dotCode (dotcode.c:1777) > ==10878== by 0x80C45A9: Rf_eval (eval.c:444) > ==10878== by 0x80C600A: do_set (eval.c:1340) > ==10878== by 0x80C44A2: Rf_eval (eval.c:427) > ==10878== by 0x80C6097: do_begin (eval.c:1104) > ==10878== by 0x80C44A2: Rf_eval (eval.c:427) > ==10878== by 0x80C61AD: do_repeat (eval.c:1066) > ==10878== by 0x80C44A2: Rf_eval (eval.c:427) > ==10878== by 0x80C6097: do_begin (eval.c:1104) > > > They also both go on to crash so hard that Valgrind crashes and says > > --9895-- INTERNAL ERROR: Valgrind received a signal 11 (SIGSEGV) - exiting > --9895-- si_code=1 Fault EIP: 0xB00313AD; Faulting address: 0x3C439F95 > --9895-- esp=0xB0653E48 > valgrind: the `impossible' happened: > Killed by fatal signal > > > Ouch.
So that looks like an nlme problem, but there's nothing obviously wrong in the code. Doug, can you spot what the problem might be? Duncan Murdoch > > -thomas > > > >> # For corExp: >> >> set.seed(26) >> for(i in 1:length(CO2$conc)){ >> CO2$conc[i]<-(CO2$conc[i]+rnorm(1)) >> } >> >> fm1CO2.lis <- nlsList(SSasympOff, CO2) >> fm1CO2.nlme <- nlme(fm1CO2.lis, control = list(tolerance = 1e-2)) >> fm2CO2.nlme <- update(fm1CO2.nlme, random = Asym + lrc ~ 1) >> CO2.nlme.var <- update(fm2CO2.nlme, >> fixed = list(Asym ~ Type * Treatment, lrc + c0 ~ 1), >> start = c(32.412, 0, 0, 0, -4.5603, 49.344), >> weights=varConstPower(fixed=list(const=0.1, power=1)), verbose=T) >> >> CO2.nlme.exp<-update(CO2.nlme.var, >> correlation=corExp(form=~as.numeric(conc)|Plant,nugget=F), data=CO2) >> >> # For corGaus: >> >> set.seed(4) >> for(i in 1:length(CO2$conc)){ >> CO2$conc[i]<-(CO2$conc[i]+rnorm(1)) >> } >> >> fm1CO2.lis <- nlsList(SSasympOff, CO2) >> fm1CO2.nlme <- nlme(fm1CO2.lis, control = list(tolerance = 1e-2)) >> fm2CO2.nlme <- update(fm1CO2.nlme, random = Asym + lrc ~ 1) >> CO2.nlme.var <- update(fm2CO2.nlme, >> fixed = list(Asym ~ Type * Treatment, lrc + c0 ~ 1), >> start = c(32.412, 0, 0, 0, -4.5603, 49.344), >> weights=varConstPower(fixed=list(const=0.1, power=1)), verbose=T) >> >> CO2.nlme.gauss<-update(CO2.nlme.var, >> correlation=corGaus(form=~as.numeric(conc)|Plant,nugget=F), data=CO2) >> >> >> On Thu, 20 Jul 2006, Duncan Murdoch wrote: >> >>> On 7/18/2006 2:28 PM, [EMAIL PROTECTED] wrote: >>>> I am having trouble fitting correlation structures within nlme. I would >>>> like to fit corCAR1, corGaus and corExp correlation structures to my data. >>>> I either get the error "step halving reduced below minimum in pnls step" >>>> or alternatively R crashes. >>>> >>>> My dataset is similar to the CO2 example in the nlme package. The one >>>> major difference is that in my case the 'conc' steps are not the same for >>>> each 'Plant'. I have replicated the problem using the CO2 data in nlme >>>> (based off of the Ch08.R script). >>>> >>>> This works (when 'conc' is the same for each 'Plant': >>>> >>>> (fm1CO2.lis <- nlsList(SSasympOff, CO2)) >>>> (fm1CO2.nlme <- nlme(fm1CO2.lis, control = list(tolerance = 1e-2))) >>>> (fm2CO2.nlme <- update(fm1CO2.nlme, random = Asym + lrc ~ 1)) >>>> CO2.nlme.var <- update(fm2CO2.nlme, >>>> fixed = list(Asym ~ Type * Treatment, lrc + c0 ~ 1), >>>> start = c(32.412, 0, 0, 0, -4.5603, 49.344), >>>> weights=varConstPower(fixed=list(const=0.1, power=1)), verbose=T) >>>> >>>> CO2.nlme.CAR<-update(CO2.nlme.var, corr=corCAR1()) >>>> >>>> CO2.nlme.gauss<-update(CO2.nlme.var, >>>> correlation=corGaus(form=~as.numeric(conc)|Plant,nugget=F), data=CO2) >>>> >>>> CO2.nlme.exp<-update(CO2.nlme.var, >>>> correlation=corExp(form=~as.numeric(conc)|Plant,nugget=F), data=CO2) But, >>>> if i change each of the 'conc' numbers slightly so that they are no longer >>>> identical between subjects i can only get the corCAR1 correlation to work >>>> while R crashes for both corExp and corGaus: >>>> >>>> for(i in 1:length(CO2$conc)){ >>>> CO2$conc[i]<-(CO2$conc[i]+rnorm(1)) >>>> } >>>> >>>> (fm1CO2.lis <- nlsList(SSasympOff, CO2)) >>>> (fm1CO2.nlme <- nlme(fm1CO2.lis, control = list(tolerance = 1e-2))) >>>> (fm2CO2.nlme <- update(fm1CO2.nlme, random = Asym + lrc ~ 1)) >>>> CO2.nlme.var <- update(fm2CO2.nlme, >>>> fixed = list(Asym ~ Type * Treatment, lrc + c0 ~ 1), >>>> start = c(32.412, 0, 0, 0, -4.5603, 49.344), >>>> weights=varConstPower(fixed=list(const=0.1, power=1)), verbose=T) >>>> >>>> CO2.nlme.CAR<-update(CO2.nlme.var, corr=corCAR1()) >>>> >>>> CO2.nlme.gauss<-update(CO2.nlme.var, >>>> correlation=corGaus(form=~as.numeric(conc)|Plant,nugget=F), data=CO2) >>>> >>>> CO2.nlme.exp<-update(CO2.nlme.var, >>>> correlation=corExp(form=~as.numeric(conc)|Plant,nugget=F), data=CO2) I >>>> have read Pinheiro & Bates (2000) and i think that it should be possible >>>> to fit these correlation structures to my data, but maybe i am mistaken. >>>> >>>> I am running R 2.3.1 and have recently updated all packages. >>> I reproduced this once in R-patched, but since then have been unable to do >>> so. I can reproduce it reliably with "set.seed(1)" at the start in R 2.3.1. >>> So it looks to me as though we've probably done something to make the error >>> less likely, but not completely fixed it. >>> >>> If you can find a script (including set.seed() to some value) that reliably >>> causes a crash in R-patched, could you let me know? >>> >>> You can get R-patched from CRAN in the bin/windows/base directory. >>> >>> Duncan Murdoch >>> >> ______________________________________________ >> R-help@stat.math.ethz.ch mailing list >> https://stat.ethz.ch/mailman/listinfo/r-help >> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html >> and provide commented, minimal, self-contained, reproducible code. >> > > Thomas Lumley Assoc. Professor, Biostatistics > [EMAIL PROTECTED] University of Washington, Seattle ______________________________________________ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.