
thanks Ben for your reply. You are absolutely right. The data is "funky" because it is data from just one mesh point (of many thousands) of a brain surface model, where the fit did not converge. So seeing other postings of converge problems I thought it might be interesting to see why one opt function converges and the other not. Probably the R-mailing list does not let through attachments, I don't know. Sorry for the inconvience.

As I would like to invesigate linear trends, do you think I could "reduce" the problem by considering the factor as ordered (what from the experiment would make sense) or even fitting a slope and applying a corAR1() correlation structure ? This should reduce the amount of parameters,shouldn't it ?

Thanks,

Stephan




----- Mensaje original -----
De: Ben Bolker <bolker@zoology.ufl.edu>
Fecha: Martes, Abril 29, 2008 17:29
Asunto: Re: [R-sig-ME] nlme optim control option
A: Stephan Moratti <moratti@med.ucm.es>, R Mixed Models <r-sig-mixed-models@r-project.org>

>    I may be missing something, but I took a quick look 
> at the
> data etc..  Your RData file didn't come through as an 
> attachment 
> (perhaps you could post it somewhere); I cut-and-pasted the
> values shown here, so I have a little bit of rounding.  
> Perhaps among 
> other things corSymm() is too general/too many parameters for 
> your data?
> You're trying to estimate 8 fixed effects, 2 variance parameters,
> and n(n-1)/2 = 28 correlations from 96 data points ... on the other
> hand, the data don't look particularly funky:
> 
> 
> Stephan Moratti wrote:
> > 
> > Thank you very much for your fast response. I will first 
> describe the 
> > data and then the models applied. Further, I will attach an 
> .RData 
> > object with all the data, so you can directly use it if you like.
> > 
> > The data are amplitude values obtained from 
> magnetoencephalografic data 
> > (like EEG). There are 8 conditions and 12 subjects (I know, 
> the number 
> > of conditions is quite high with respect to number of 
> subjects). The 
> > basic data matrix looks like this (called "mat" in the RData 
> object):> 
> >         
> V1       V2       V3       V4       V5       V6       V7       V8
> > 1  2.851582 3.314363 3.279070 2.394195 2.914685 1.876838 
> 2.180715 3.203826
> > 2  4.521648 6.171458 4.124039 5.010868 4.923190 5.221255 
> 4.745925 5.477719
> > 3  3.957449 4.743048 3.015077 2.890330 3.428396 3.078638 
> 3.015924 2.921296
> > 4  5.139640 5.272191 4.067972 3.151300 5.739077 3.279299 
> 3.524670 6.110507
> > 5  7.928020 8.344273 7.525576 6.835657 8.995732 8.871082 
> 7.821745 8.346335
> > 6  4.900010 4.508470 2.532844 4.197058 4.701989 2.878392 
> 4.939342 2.913649
> > 7  4.378253 7.150406 5.205363 4.872032 5.460943 5.901771 
> 4.986221 4.985459
> > 8  5.697526 4.945624 3.855982 4.460146 4.661727 5.434473 
> 6.283447 5.795619
> > 9  2.772407 2.823246 2.693566 3.019434 2.303198 3.340165 
> 3.624334 2.974910
> > 10 5.063280 6.384693 6.015056 7.935250 5.564225 7.541611 
> 7.154139 8.031693
> > 11 3.202468 1.877924 1.716837 2.286963 2.152085 3.096542 
> 3.026434 4.369565
> > 12 4.684218 5.391280 5.360308 4.675992 4.500356 3.507834 
> 3.795734 7.622723
> > 
> > The same data matrix is transformed to a data frame ("datamat" 
> in the 
> > RData object):
> >     values ind subject condition
> > 1  2.851582  V1       
> 1         1
> > 2  4.521648  V1       
> 2         1
> > 3  3.957449  V1       
> 3         1
> > 4  5.139640  V1       
> 4         1
> > 5  7.928020  V1       
> 5         1
> > 6  4.900010  V1       
> 6         1
> > 7  4.378253  V1       
> 7         1
> > 8  5.697526  V1       
> 8         1
> > 9  2.772407  V1       
> 9         1
> > 10 5.063280  V1      
> 10         1
> > 11 3.202468  V1      
> 11         1
> > 12 4.684218  V1      
> 12         1
> > 13 3.314363  V2       
> 1         2
> > 14 6.171458  V2       
> 2         2
> > 15 4.743048  V2       
> 3         2
> > 16 5.272191  V2       
> 4         2
> > 17 8.344273  V2       
> 5         2
> > 18 4.508470  V2       
> 6         2
> > 19 7.150406  V2       
> 7         2
> > 20 4.945624  V2       
> 8         2
> > 21 2.823246  V2       
> 9         2
> > 22 6.384693  V2      
> 10         2
> > 23 1.877924  V2      
> 11         2
> > 24 5.391280  V2      
> 12         2
> > .... up to 8 conditions
> > 
> > **********************
> > * Applying optim resutls in: *
> > **********************
> > 
> > fm.optim<-
> lme(values~condition,data=datamat,random=~1|subject,correlation=corSymm(),control=list(msMaxIter=100,opt 
> > = "optim",msVerbose=TRUE))
> > initial  value 209.185062
> > iter  10 value 186.050690
> > iter  20 value 184.942325
> > iter  30 value 184.724836
> > iter  40 value 184.121038
> > iter  50 value 183.918251
> > final  value 183.903769
> > converged
> > 
> > summary(fm.optim)
> > 
> > Linear mixed-effects model fit by REML
> > Data: datamat
> >       
> AIC      BIC    logLik
> >  299.5351 393.6739 -111.7675
> > 
> > Random effects:
> > Formula: ~1 | subject
> >        (Intercept) Residual
> > StdDev:    1.110893 1.176853
> > 
> > Correlation Structure: General
> > Formula: ~1 | subject
> > Parameter estimate(s):
> > Correlation:
> >  1      
> 2      3      
> 4      5      
> 6      7    2 
> > -
> 0.056                                         3 -0.432  
> > 
> 0.742                                  4 -0.213  0.545  
> > 
> 0.586                           5  0.729  0.593  0.281  
> > 
> 0.235                    6  0.101  0.474  0.412  0.657  
> > 
> 0.430             7  0.625  0.054 -0.212  0.440  0.509  0.635      8 
> > -0.087  0.261  0.557  0.417  0.257  
> 0.354  0.005
> > Fixed effects: values ~ condition
> >                Value Std.Error DF   t-value p-value
> > (Intercept)  4.591375 0.4671782 77  9.827888  0.0000
> > condition2   0.485873 0.4936256 77  
> 0.984295  0.3281
> > condition3  -0.475401 0.5748486 77 -0.827002  0.4108
> > condition4  -0.280606 0.5290556 77 -0.530391  0.5974
> > condition5   0.020759 0.2498983 77  
> 0.083068  0.9340
> > condition6  -0.089050 0.4556203 77 -0.195448  0.8456
> > condition7   0.000178 0.2941852 77  
> 0.000604  0.9995
> > condition8   0.638067 0.5009446 77  
> 1.273727  0.2066
> > Correlation:
> >           
> (Intr) cndtn2 cndtn3 cndtn4 cndtn5 cndtn6 cndtn7
> > condition2 -
> 0.528                                         condition3 
> > -0.615  
> 0.907                                  condition4 -0.566  0.801  
> > 
> 0.846                           condition5 -0.267  0.860  0.790  
> > 
> 0.627                    condition6 -0.488  0.734  0.768  0.847  
> > 
> 0.608             condition7 -0.315  0.385  0.406  0.762  0.243  
> > 0.783      condition8 -0.536  
> 0.655  0.832  0.748  0.566  0.678  0.366
> > 
> > Standardized Within-Group Residuals:
> >       
> Min         Q1        Med         Q3        Max
> > -1.9205059 -0.7430723 -0.1023267  0.7054844  2.1579685
> > 
> > Number of Observations: 96
> > Number of Groups: 12
> > 
> > ***********************
> > * Applying nlminb results in: *
> > ***********************
> > fm.nlminb<-
> lme(values~condition,data=datamat,random=~1|subject,correlation=corSymm(),control=list(msMaxIter=100,opt 
> > = "nlminb",msVerbose=TRUE))
> > 
> > ....
> > 98:     183.73454: 0.127036 0.0579634 
> 0.506547 -1.36663 0.306432 
> > -0.838546 -0.395783 -1.18434 -2.01708 -1.91615 0.0747060 -
> 0.199496 
> > -0.740317 -0.398701 -0.835580 0.623620 -0.979961 -0.232816 -
> 0.0991876 
> > -1.63724  2.08835 -12.9687 0.147038 -0.361111 -0.964921 -
> 0.271664 
> > -0.512183 -0.409421 -2.67278
> > 99:     183.73454: 0.127040 0.0579607 
> 0.506550 -1.36663 0.306427 
> > -0.838537 -0.395789 -1.18434 -2.01709 -1.91611 0.0747092 -
> 0.199503 
> > -0.740320 -0.398707 -0.835577 0.623649 -0.979960 -0.232836 -
> 0.0991891 
> > -1.63723  2.08842 -12.9657 0.147034 -0.361111 -0.964925 -
> 0.271658 
> > -0.512177 -0.409425 -2.67234
> > 100:     183.73454: 0.127046 0.0579570 
> 0.506540 -1.36663 0.306412 
> > -0.838550 -0.395794 -1.18434 -2.01711 -1.91609 0.0747311 -
> 0.199514 
> > -0.740328 -0.398712 -0.835571 0.623651 -0.979969 -0.232830 -
> 0.0992005 
> > -1.63723  2.08843 -12.9642 0.147026 -0.361116 -0.964933 -
> 0.271660 
> > -0.512178 -0.409420 -2.67211
> > 100:     183.73454: 0.127046 0.0579570 
> 0.506540 -1.36663 0.306412 
> > -0.838550 -0.395794 -1.18434 -2.01711 -1.91609 0.0747311 -
> 0.199514 
> > -0.740328 -0.398712 -0.835571 0.623651 -0.979969 -0.232830 -
> 0.0992005 
> > -1.63723  2.08843 -12.9642 0.147026 -0.361116 -0.964933 -
> 0.271660 
> > -0.512178 -0.409420 -2.67211
> > Error en lme.formula(values ~ condition, data = datamat, 
> random = ~1 |  :
> >  nlminb problem, convergence error code = 1
> >  message = iteration limit reached without convergence (9)
> > 
> > using more iterations:
> > fm.nlminb2<-
> lme(values~condition,data=datamat,random=~1|subject,correlation=corSymm(),control=list(msMaxIter=200,opt 
> > ="nlminb",msVerbose=TRUE))
> > 
> > ....
> > 98:     183.73454: 0.127036 0.0579634 
> 0.506547 -1.36663 0.306432 
> > -0.838546 -0.395783 -1.18434 -2.01708 -1.91615 0.0747060 -
> 0.199496 
> > -0.740317 -0.398701 -0.835580 0.623620 -0.979961 -0.232816 -
> 0.0991876 
> > -1.63724  2.08835 -12.9687 0.147038 -0.361111 -0.964921 -
> 0.271664 
> > -0.512183 -0.409421 -2.67278
> > 99:     183.73454: 0.127040 0.0579607 
> 0.506550 -1.36663 0.306427 
> > -0.838537 -0.395789 -1.18434 -2.01709 -1.91611 0.0747092 -
> 0.199503 
> > -0.740320 -0.398707 -0.835577 0.623649 -0.979960 -0.232836 -
> 0.0991891 
> > -1.63723  2.08842 -12.9657 0.147034 -0.361111 -0.964925 -
> 0.271658 
> > -0.512177 -0.409425 -2.67234
> > 100:     183.73454: 0.127046 0.0579570 
> 0.506540 -1.36663 0.306412 
> > -0.838550 -0.395794 -1.18434 -2.01711 -1.91609 0.0747311 -
> 0.199514 
> > -0.740328 -0.398712 -0.835571 0.623651 -0.979969 -0.232830 -
> 0.0992005 
> > -1.63723  2.08843 -12.9642 0.147026 -0.361116 -0.964933 -
> 0.271660 
> > -0.512178 -0.409420 -2.67211
> > 101:     183.73454: 0.127045 0.0579574 
> 0.506540 -1.36663 0.306412 
> > -0.838549 -0.395794 -1.18434 -2.01711 -1.91609 0.0747313 -
> 0.199515 
> > -0.740328 -0.398712 -0.835570 0.623651 -0.979969 -0.232834 -
> 0.0991995 
> > -1.63723  2.08843 -12.9642 0.147026 -0.361116 -0.964934 -
> 0.271660 
> > -0.512178 -0.409419 -2.67211
> > 102:     183.73454: 0.127045 0.0579574 
> 0.506540 -1.36663 0.306412 
> > -0.838549 -0.395794 -1.18434 -2.01711 -1.91609 0.0747313 -
> 0.199515 
> > -0.740328 -0.398712 -0.835570 0.623651 -0.979969 -0.232834 -
> 0.0991995 
> > -1.63723  2.08843 -12.9642 0.147026 -0.361116 -0.964934 -
> 0.271660 
> > -0.512178 -0.409419 -2.67211
> > Error en lme.formula(values ~ condition, data = datamat, 
> random = ~1 |  :
> >  nlminb problem, convergence error code = 1
> >  message = false convergence (8)
> > 
> > However, if I log transform the data, both optimization 
> methods converge 
> > using a maxIter of 200. All objects could be found in the 
> Rdata object. 
> > Not to blow up the mail, I did not print the summary of the 
> two fitted 
> > models as they are in the RData object. The two methods 
> converge at 
> > similar values using the log transformed data. I don't know 
> why "nlminb" 
> > does not converge with the not transformed data.
> > 
> > Best,
> > 
> > Stephan
> > 
> > 
> > 
> > 
> > Douglas Bates escribió:
> >> On 4/29/08, Doran, Harold <HDoran@air.org> wrote:
> >>  
> >>> I actually trust nlminb more than I trust optim. In some of 
> the 
> >>> functions I have written in the MiscPsycho package (e.g., 
> >>> irt.ability) we get multimodal distributions for some of the 
> item 
> >>> response theory models.
> >>>
> >>>  When I use optim for optimization, you can get very 
> funny results 
> >>> with difficult to maximize distributions. nlminb, however, 
> gave 
> >>> reasonable results and the maximum was always confirmed to 
> be 
> >>> "correct" via a visual examination of the likelihoods.
> >>>
> >>>  In other words, I can plot the likelihood function and 
> get a sense 
> >>> of where the max is. optim would more often than not give 
> results 
> >>> that are very (very) far away from the max with multimodal 
> >>> distributions whereas nlminb always gave back the max that 
> seemed 
> >>> consistent with the visual plot of the likelihood.
> >>>
> >>>  With all of this said, you should switch to lme4.
> >>>     
> >>
> >> As Harold indicated, switching the optimizer used in the nlme 
> package>> from optim to nlminb was intentional.  I had seen 
> cases of spurious
> >> convergence for optim and found it easier to introduce the nlminb
> >> optimizer.
> >>
> >> Could you give more details of the cases where functions from 
> the nlme
> >> package using nlminb are not converging?  Are you using 
> the lme
> >> function or the nlme function?  Does your model have 
> random effects
> >> and fixed effects only or does it also use parameterized 
> correlation>> structures?
> >>
> >> As Harold indicated, the lme4 package is now the preferred 
> way to fit
> >> linear or nonlinear mixed-effects models that have only fixed 
> effects>> and random effects.  If you can show us one of 
> your calls to lme or
> >> nlme we can tell you how to translate it to a call to lmer or 
> nlmer.>>
> >>  
> >>>  > -----Original Message-----
> >>>  > From: r-sig-mixed-models-bounces@r-project.org
> >>>  > [mailto:r-sig-mixed-models-bounces@r-project.org] On 
> Behalf>>>  > Of Stephan Moratti
> >>>  > Sent: Tuesday, April 29, 2008 7:05 AM
> >>>  > To: r-sig-mixed-models@r-project.org
> >>>  > Subject: [R-sig-ME] nlme optim control option
> >>>  >
> >>>  > Hello,
> >>>  >
> >>>  > I am using the nlme package for quite a while, but 
> recently I
> >>>  > have discovered the "optim" option in the lmeControl 
> function>>>  > and I have read that until the R verison 
> 2.2.0 the
> >>>  > optimization function is "nlminb" and before it was 
> "optim">>>  > by default. As I am a "user" of R I am not 
> into the
> >>>  > optimization functions. However, I have a case now 
> where with
> >>>  > optim = "nlminb" the model does not converge, 
> whereas with
> >>>  > optim = "optim" the model converges.
> >>>  >
> >>>  > Can somebody explain to me in easy words what the 
> difference>>>  > is ? If nlminb does not converge, but 
> optim does, is the
> >>>  > result less trustable ?
> >>>  >
> >>>  > Thanks,
> >>>  >
> >>>  > Stephan
> >>>  >
> >>>  >
> >>>  > --
> >>>  > *Stephan Moratti, PhD/
> >>>  > /**/see also: http://web.mac.com/smoratti/ /*Centro de
> >>>  > Tecnología Biomédica CBT, Universidad Politécnica de 
> Madrid,>>>  >
> >>>  > en la actualidad (currently at) en el
> >>>  > Centro de Magnetoencefalografía Dr. Perez Modrego,
> >>>  > Universidad Complutense de Madrid, Faculdad de Medicina,
> >>>  > Pabellón 8, Avda. Complutense, s/n, 28040 Madrid, Spain,
> >>>  >
> >>>  > email: moratti@gbt.tfo.upm.es 
> <mailto:moratti@gbt.tfo.upm.es>>>>  
> >          moratti@med.ucm.es
> >>>  > Tel.:    +34 91 394 2292
> >>>  > Fax.:   +34 91 394 2294
> >>>  > */
> >>>  > /*
> >>>  >
> >>>  > _______________________________________________
> >>>  > R-sig-mixed-models@r-project.org mailing list
> >>>  > https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> >>>  >
> >>>
> >>>  _______________________________________________
> >>>  R-sig-mixed-models@r-project.org mailing list
> >>>  https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> >>>
> >>>     
> >>
> >>   
> > 
> > 
> > 
> > ---------------------------------------------------------------
> ---------
> > 
> > _______________________________________________
> > R-sig-mixed-models@r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> 
> 
>

	[[alternative HTML version deleted]]


