[R-es] Ajuste no lineal con nls no converge

Argel Gastélum Arellánez argel.gastelum en gmail.com
Mar Ene 24 01:26:09 CET 2012


     Hola Olivier, muchas gracias por tus comentarios. Pongo algunas 
notas más abajo.

El 23/01/12 16:43, Olivier Nuñez escribió:
> Argel,
>
> Kiu = 1.2 x 10^64 es un valor bastante grande para una constante de 
> inhibición!

    Efectivamente, Kiu = 1.2 x 10^64 es demasiado grande para una 
constante de inhibición. El problema no es de interpretación del 
resultado, es más del por qué "nls" no converge hacia los mismos valores 
que ellos reportan para el mismo set de datos crudos, es más, sólo no 
converge. Es precisamente por lo que los autores del artículo ajustan a 
los mismos datos el modelo que considera sólo la inhibición competitiva.

> Y sospecho que el algoritmo de estos autores tampoco converge con 
> estos datos.

     Los autores mencionan que en SAS encuentran convergencia con este 
modelo mixto hacia los valores que mencioné anteriormente. Yo realmente 
desconozco SAS y no tengo acceso a él como para comprobarlo. Pensé que 
con R y "nls" debería llegar a converger también puesto que con 
algorithm = "default" se trataría también del método Gauss-Newton que 
ellos usaron en SAS.

> De hecho si pones trace=TRUE, verás que este parámetro se va al 
> infinito y el problema parece más bien de identificación del modelo:
>
     Ya lo había hecho antes de consultar a la lista y había notado que 
Kiu se va incrementando hacia valores muy grandes y sin mucho sentido 
(el comportamiento esperado, ya que el modelo mixto no se ajusta bien a 
los datos, sino que se ajusta mejor a un modelo con sólo inhibición 
competitiva), pero como los autores mencionan convergencia creí que en R 
también pasaría así. Lo que me interesa por el momento es más bien poder 
obtener los mismos valores ajustados que ellos reportan para los 
parámetros, sabiendo de antemano que no son los que mejor se ajustan 
(explico esto con un poco más de detalle más abajo).

> > AJUSTE.default <- nls(t ~ 
> ((Km+Km*P0/Kic+Km*S0/Kic)*log((S0-(Pt-P0))/S0)+(-Km/Kic+1+P0/Kiu+S0/Kiu)*((S0-(Pt-P0))-S0)+(-1/(2*Kiu))*((S0-(Pt-P0))^2-S0^2))*(-1/(E*k3)), 
> DATOS, start = list (Km = 50, Kic = 6, Kiu = 6, k3 = 0.02), algorithm 
> = "default",trace=TRUE)
> 9689.914 :  50.00  6.00  6.00  0.02
> 6670.874 :  34.93095937  5.79629581 14.16884101  0.01610926
> 1896.282 :  37.29205443  5.77977412 32.30063332  0.01671887
> 840.2903 :  37.3846129  5.7802131 88.8218065  0.0167428
> 536.3691 :   37.38476020   5.78018175 361.40233750   0.01674285
> 463.0871 :  3.738474e+01 5.780188e+00 3.765039e+03 1.674285e-02
> 455.831 :  3.738453e+01 5.780068e+00 3.377163e+05 1.674283e-02
> 455.8272 :  3.738363e+01 5.780040e+00 3.325798e+08 1.674261e-02
> Error en nls(t ~ ((Km + Km * P0/Kic + Km * S0/Kic) * log((S0 - (Pt - 
> P0))/S0) +  :
>   singular gradient
>
> Un saludo. Olivier

     Comento un poco más sobre el artículo. El objetivo en general es 
ajustar los siguientes modelos de Michaelis-Menten: Sin inhibición, con 
inhibición competitiva, con inhibición no competitiva, con inhibición 
acompetitiva y con inhibición mixta (todos en sus formas integradas), y 
de acuerdo con los resultados de los valores de los parámetros, del 
número de parámetros y de las sumas de cuadrados decidir cuál es el que 
mejor se ajusta a los datos.

     En R, "nls" ha reproducido con exactitud los parámetros estimados 
por los autores para los modelos sin inhibición y con inhibición 
competitiva, pero no para los que consideran inhibición acompetitiva e 
inhibición mixta (donde "nls" reporta no convergencia y los autores 
reportan convergencia con valores de Kiu muy grandes, del orden de 10^38 
a 10^64). Los autores mencionan que ellos tampoco encontraron 
convergencia con el modelo de inhibición no competitiva. Sin embargo, al 
usar en R y "nls" las opciones warnOnly = TRUE, "trrosace" y "port" 
observo que aunque no reporta convergencia los valores de los parámetros 
ajustados son muy similares a lo reportado por los autores, excepto para 
las constantes de inhibición que resultan demasiado grandes, donde "nls" 
reporta no convergencia y valores finales del orden de 10^5 para Kiu.

     Espero haber aclarado un poco la falta de información preliminar.

     De nuevo muchas gracias por tus comentarios.

     Saludos.

--
     Argel.



Más información sobre la lista de distribución R-help-es