FANDOM


INFORME SOBRE AJUSTE NO LINEALEditar

Autores: Carlos Mohor, Pablo Rios, Alfredo Valenzuela

Curso: Métodos Numéricos

Profesor: Stefan Berres

Universidad Católica de Temuco

Temuco, 27 de Noviembre del 2012

Resumen EjecutivoEditar

El presente informe expone el tema de análisis de modelos no-lineales. 

El objetivo es que, dado un modelo no-lineal,  con variables desconocidas, podamos encontrar los parámetros que ajusten el modelo de mejor manera a los datos observados.

Una vez que podamos entender bien la lógica, la idea es poder formular el algoritmo para implementarlo en código y analizar los resultados.

Para esto hemos estructurado el informe de la siguiente manera:

I.                    Introducción

II.                  Explicación paso a paso de cómo buscar los parámetros, linealizar el modelo y obtener las ecuaciones en base a un modelo de función logística trabajado en clases.

III.                Formulación de un algoritmo general en MATLAB

IV.                Implementación paso a paso de un algoritmo para un modelo no-lineal distinto (variante).

V.                  Implementación paso a paso del algoritmo pero en un modelo nuevo (no visto en clases), previa introducción.

I. IntroducciónEditar

¿Qué es un modelo matemático?


• Conjunto de ecuaciones que relacionan las variables de interés del proceso y representan adecuadamente su comportamiento

• Siempre son aproximaciones de la realidad

• Distintos modelos para distintos objetivos y tipos de procesos

• Compromiso entre facilidad de uso y exactitud

"NO LINEALIDAD" 

Los sistemas no lineales representan sistemas cuyo comportamiento no es expresable como la suma de los comportamientos de sus descriptores. Más formalmente, un sistema físico, matemático o de otro tipo es no lineal cuando las ecuaciones de movimiento, evolución o comportamiento que regulan su comportamiento son no lineales. En particular, el comportamiento de sistemas no lineales no está sujeto al principio de superposición, como lo es un sistema lineal.

La linealidad de un sistema permite a los investigadores hacer ciertas suposiciones matemáticas y aproximaciones, permitiendo un cálculo más sencillo de los resultados. Ya que los sistemas no lineales no son iguales a la suma de sus partes, usualmente son difíciles (o imposibles) de modelar, y sus comportamientos con respecto a una variable dada (por ejemplo, el tiempo) es extremadamente difícil de predecir.

Algunos sistemas no lineales tienen soluciones exactas o integrables, mientras que otros tienen comportamiento caótico, por lo tanto no se pueden reducir a una forma simple ni se pueden resolver. Un ejemplo de comportamiento caótico son las olas gigantes. Aunque algunos sistemas no lineales y ecuaciones de interés general han sido extensamente estudiados, la vasta mayoría son pobremente comprendidos.

II. Explicación paso a paso el desarrollo del modelo de función logísticaEditar

OBJETIVO: Entender la lógica del desarrollo

dado el modelo:

f(x)=\frac{a}{b+c\exp(-dx)}

y los datos

\begin{array}{|c|c|c|c|c|c|c|}
      \hline
      X_i & 1 & 2 & 3 & 4 & 5 & 6\\
      \hline
      f_i & 1 & 4 & 15 & 30 & 50 & 55\\
      \hline
      \end{array}

formular el algoritmo que minimiza

S= \sum_{i=1}^n [f_i-f(x_i)]^2

buscar los parametros a, b, c, d que minimiza S

Criterio de extreme(minimo)

\frac{ds}{da} = 0

\frac{ds}{db} = 0

\frac{ds}{dc} = 0

\frac{ds}{dd} = 0

Tenemos 4 parametros (como variables desconocidas) y 6 ecuaciones(que corresponden a (x_i,f_i), i=1, ... , 6))

Linealizar con respecto a \begin{array}{llll}\Delta a, \Delta b, \Delta c, \Delta d \end{array}

Buscar y calcular los parametros a^*, b^*, c^*, d^* , tal que

S(a^*, b^*, c^*, d^*)\min S(a, b, c, d)

Parametros a, b, c, d para que la curva (x,f(x)) se ajuste de mejor manera a los datos.

Los parametros deben coincidir con la escala, se ajustan a la cuerva(podrian tener la misma forma, pero en otra escala).

Paso1:Editar

"Transformar eñ problema de minimizacion a una tarea de calcular ceros"

Queremos minimizar S(a, b, c, d)


\begin{array}{llll}
\left . 
       \begin{matrix} 
          F_1(a,b,c,d):=\frac{\partial S(a,b,c,d)}{\partial a}=0\\
          F_2(a,b,c,d):=\frac{\partial S(a,b,c,d)}{\partial b}=0\\
          F_3(a,b,c,d):=\frac{\partial S(a,b,c,d)}{\partial c}=0\\
          F_4(a,b,c,d):=\frac{\partial S(a,b,c,d)}{\partial d}=0\\
       \end{matrix} 
    \right \}=F=\Delta S(a,b,c,d)\end{array}

Definir: Todos los componentes de la funcion F dependen de los parametros(dependencia).

F = \left ( 
      \begin{matrix} 
         F_1 \\
         F_2 \\
         F_3 \\
         F_4
      \end{matrix}
   \right ) , F = F(a,b,c,d)

\begin{array}{llll}F =\Delta S; \end{array} (F corresponde a la gradiente de la funcion S).

Para encontrar el minimo de S(a,b,c,d) necesitamos encontrar el cero de F:

F(a,b,c,d) =\left ( 
      \begin{matrix} 
         0 \\
         0 \\
         0 \\
         0
      \end{matrix}
   \right )

Calculando los componentes de F:

F_1(a,b,c,d) =\frac{\partial S(a,b,c,d)}{\partial a}

  = \frac{\partial }{\partial a} \sum_{i=1}^n (f(x_i)-f_i)^2)

  =  \sum_{i=1}^n \frac{\partial }{\partial a}(f(x_i)-f_i)^2)

  =  \sum_{i=1}^n  2(f(x_i)-f_i)  \frac{\partial }{\partial a} f(x_i)

... lo mismo para los otros F_2, F_3, F_4 ...

Por lo tanto;

F_1(a,b,c,d) = \left [  \sum_{i=1}^n  2(f(x_i)-f_i)  \frac{\partial }{\partial a} f(x_i) \right ]

F_2(a,b,c,d) = \left [  \sum_{i=1}^n  2(f(x_i)-f_i)  \frac{\partial }{\partial a} f(x_i) \right ]

F_3(a,b,c,d) = \left [  \sum_{i=1}^n  2(f(x_i)-f_i)  \frac{\partial }{\partial a} f(x_i) \right ]

F_4(a,b,c,d) = \left [  \sum_{i=1}^n  2(f(x_i)-f_i)  \frac{\partial }{\partial a} f(x_i) \right ]

Paso2:Editar

"Linealizar el modelo no lineal"

Dificultad: -Entender y aplicar rotacion vectorial. -Entender la secuancia de varios pasos de transformacion del problema

Linealizar el modelo no lineal es decir, linealizar:

F(a,b,c,d)

Para aproximar F(P)=0 , con P=(a,b,c,d) ,

Introduciremos P_i=(a_i,b_i,c_i,d_i) , i=0,1,2,...

Formular el metodo de newton, con la especificacion de la funcion:

F(a,b,c,d)


F(P_{i+1})=F(p_i)+\left ( 
      \begin{matrix} 
         \frac{\partial F_1}{\partial a} & \frac{\partial F_1}{\partial b} & \frac{\partial F_1}{\partial c} & \frac{\partial F_1}{\partial d}\\
         \frac{\partial F_2}{\partial a} & \frac{\partial F_2}{\partial b} & \frac{\partial F_2}{\partial c} & \frac{\partial F_2}{\partial d}\\
         \frac{\partial F_3}{\partial a} & \frac{\partial F_3}{\partial b} & \frac{\partial F_3}{\partial c} & \frac{\partial F_3}{\partial d}\\
         \frac{\partial F_4}{\partial a} & \frac{\partial F_4}{\partial b} & 
\frac{\partial F_4}{\partial c} & \frac{\partial F_4}{\partial d}
      \end{matrix}
   \right ) \Delta P = 0

Iteracion con respecto a i \to i+1

(1) g(P_i)\Delta P = -F(p_i) \to con esto podremos implementar el algoritmo

(2) P_{i+1}=P_i+\Delta P

donde \Delta P=P_{i+1}-P_i=\left ( 
      \begin{matrix} 
         a_{i+1}- a_i \\
         b_{i+1}- b_i \\
         c_{i+1}- c_i \\
         d_{i+1}- d_i
      \end{matrix}
   \right )=\left ( 
      \begin{matrix} 
         \Delta a \\
         \Delta b \\
         \Delta c \\
         \Delta d
      \end{matrix}
   \right )

-Especificar la matriz g(P_i) , en termino de solamente a, b, c, d o sea:

-Especificar g(P_i) en terminos de f ,\frac{\partial f}{\partial a},\frac{\partial f}{\partial b},\frac{\partial f}{\partial c},\frac{\partial f}{\partial d}

-Calcular derivadas parciales de f .

g(P_i) , en terminos de f,a ,b ,c ,d:

Tenemos; g(p_i)=\left ( 
      \begin{matrix} 
         g_{11}(p_i) &  g_{12}(p_i) &  g_{13}(p_i) &  g_{14}(p_i)\\
         g_{21}(p_i) &  g_{22}(p_i) &  g_{23}(p_i) &  g_{24}(p_i)\\
         g_{31}(p_i) &  g_{32}(p_i) &  g_{33}(p_i) &  g_{34}(p_i)\\
         g_{41}(p_i) &  g_{42}(p_i) &  g_{43}(p_i) &  g_{44}(p_i)
      \end{matrix}
   \right )



      \begin{matrix} 
         g_{11}(p_i)\frac{\partial F_1}{\partial a} &  g_{12}(p_i)\frac{\partial F_1}{\partial b} &  g_{13}(p_i)\frac{\partial F_1}{\partial c} &  g_{14}(p_i)\frac{\partial F_1}{\partial c}\\
         g_{21}(p_i)\frac{\partial F_1}{\partial a} &  g_{22}(p_i)\frac{\partial F_1}{\partial b} &  g_{23}(p_i)\frac{\partial F_1}{\partial c} &  g_{24}(p_i)\frac{\partial F_1}{\partial c}\\
         g_{31}(p_i)\frac{\partial F_1}{\partial a} &  g_{32}(p_i)\frac{\partial F_1}{\partial b} &  g_{33}(p_i)\frac{\partial F_1}{\partial c} &  g_{34}(p_i)\frac{\partial F_1}{\partial c}\\
         g_{41}(p_i)\frac{\partial F_1}{\partial a} &  g_{42}(p_i)\frac{\partial F_1}{\partial b} &  g_{43}(p_i)\frac{\partial F_1}{\partial c} &  g_{44}(p_i)\frac{\partial F_1}{\partial c}
      \end{matrix}


\to calculamos:

(1) Para g_{11}(P_i) ;

g_{11}(P_i) = \frac{\partial F_1}{\partial a}; \to se aplica la definicion de:

F_1=f_1(a,b,c,d):= \frac{\partial S}{\partial a}(a,b,c,d)

\therefore g_{11}(P_i)= \frac{\partial F_1}{\partial a}

=\frac{\partial }{\partial a}\left [  \sum_{i=1}^n  2(f(x_i)-f_i)  \frac{\partial }{\partial a} f(x_i) \right ]

= \sum_{i=1}^n \frac{\partial }{\partial a}\left [2(f(x_i)-f_i)  \frac{\partial }{\partial a} f(x_i) \right ]

=\sum_{i=1}^n \left [2(f(x_i)-f_i)\frac{\partial }{\partial a}\left (\frac{\partial }{\partial a}f(x_i)\right )+\frac{\partial }{\partial a}f(x_i)\frac{\partial }{\partial a} 2(f(x_i)-f_i)\right ]

=\left [2(f(x_i)-f_i)\frac{\partial ^2}{\partial a^2}f(x_i) + 2\left (\frac{\partial }{\partial a}f(x_i)\right )^2\right ]

g_{11}(P_i) =\sum_{i=1}^n 2\left [f(x_i)\frac{\partial ^2}{\partial a^2}f(x_i) + \left (\frac{\partial }{\partial a}f(x_i)\right )^2\right ]


(2) Para g_{12}(P_i) ;

\frac{\partial F_1}{\partial b}=2f(x_i)\frac{\partial }{\partial b}\left (\frac{\partial }{\partial a}f(x_i)\right )+\frac{\partial }{\partial a} f(x_i)\frac{\partial }{\partial b}2f(x_i)


(3) Para g_{13}(P_i) ;

\frac{\partial F_1}{\partial c}=2f(x_i)\frac{\partial }{\partial c}\left (\frac{\partial }{\partial a}f(x_i)\right )+\frac{\partial }{\partial a} f(x_i)\frac{\partial }{\partial c}2f(x_i)


(4) Para g_{14}(P_i) ;

\frac{\partial F_1}{\partial d}=2f(x_i)\frac{\partial }{\partial d}\left (\frac{\partial }{\partial a}f(x_i)\right )+\frac{\partial }{\partial a} f(x_i)\frac{\partial }{\partial d}2f(x_i)


(5) Para g_{21}(P_i) ;

\frac{\partial F_2}{\partial a}=2f(x_i)\frac{\partial }{\partial a}\left (\frac{\partial }{\partial b}f(x_i)\right )+\frac{\partial }{\partial b} f(x_i)\frac{\partial }{\partial a}2f(x_i)


(6) Para g_{22}(P_i) ;

\frac{\partial F_2}{\partial b}=2\left [f(x_i)\frac{\partial ^2}{\partial b^2}f(x_i) + \left (\frac{\partial }{\partial b}f(x_i)\right )^2\right ]


(7) Para g_{23}(P_i) ;

\frac{\partial F_2}{\partial c}=2f(x_i)\frac{\partial }{\partial c}\left (\frac{\partial }{\partial b}f(x_i)\right )+\frac{\partial }{\partial b} f(x_i)\frac{\partial }{\partial c}2f(x_i)


(8) Para g_{24}(P_i) ;

\frac{\partial F_2}{\partial d}=2f(x_i)\frac{\partial }{\partial d}\left (\frac{\partial }{\partial b}f(x_i)\right )+\frac{\partial }{\partial b} f(x_i)\frac{\partial }{\partial d}2f(x_i)


(9) Para g_{31}(P_i) ;

\frac{\partial F_3}{\partial a}=2f(x_i)\frac{\partial }{\partial a}\left (\frac{\partial }{\partial c}f(x_i)\right )+\frac{\partial }{\partial c} f(x_i)\frac{\partial }{\partial a}2f(x_i)


(10) Para g_{32}(P_i) ;

\frac{\partial F_3}{\partial b}=2f(x_i)\frac{\partial }{\partial b}\left (\frac{\partial }{\partial c}f(x_i)\right )+\frac{\partial }{\partial c} f(x_i)\frac{\partial }{\partial b}2f(x_i)


(11) Para g_{33}(P_i) ;

\frac{\partial F_3}{\partial c}=2\left [f(x_i)\frac{\partial ^2}{\partial c^2}f(x_i) + \left (\frac{\partial }{\partial c}f(x_i)\right )^2\right ]


(12) Para g_{34}(P_i) ;

\frac{\partial F_3}{\partial d}=2f(x_i)\frac{\partial }{\partial d}\left (\frac{\partial }{\partial c}f(x_i)\right )+\frac{\partial }{\partial c} f(x_i)\frac{\partial }{\partial d}2f(x_i)


(13) Para g_{41}(P_i) ;

\frac{\partial F_4}{\partial a}=2f(x_i)\frac{\partial }{\partial a}\left (\frac{\partial }{\partial d}f(x_i)\right )+\frac{\partial }{\partial d} f(x_i)\frac{\partial }{\partial a}2f(x_i)


(14) Para g_{42}(P_i) ;

\frac{\partial F_4}{\partial b}=2f(x_i)\frac{\partial }{\partial b}\left (\frac{\partial }{\partial d}f(x_i)\right )+\frac{\partial }{\partial d} f(x_i)\frac{\partial }{\partial b}2f(x_i)


(15) Para g_{43}(P_i) ;

\frac{\partial F_4}{\partial c}=2f(x_i)\frac{\partial }{\partial c}\left (\frac{\partial }{\partial d}f(x_i)\right )+\frac{\partial }{\partial d} f(x_i)\frac{\partial }{\partial c}2f(x_i)


(16) Para g_{44}(P_i) ;

\frac{\partial F_4}{\partial c}=2\left [f(x_i)\frac{\partial ^2}{\partial d^2}f(x_i) + \left (\frac{\partial }{\partial d}f(x_i)\right )^2\right ]

III. Formulación de un algoritmo general en MATLABEditar

Seudocódigo para una implementación en MATLAB:

Basado en el modelo f(a,b;x)

Sumatoria “S” es la misma que hemos utilizado hasta ahora, con los datos Xi,Fi

%Declaracion de variables

%generacion de datos (x(i),f(i))

a=1, b=2

for i=1 hasta n %generamos los datos

x(i) = i/n;

f(i) = a*b*x(i)

end for

a=0, b=1 %escogemos los parametros iniciales de iteración

i=0; f=[1;1 ]; épsilon = 1e-7; xx=[]

%Iteracion:

While norm(F)<épsilon && i<n

%calculamos F1 (la derivada de S  con respect a A)

%f(a,b;x)=abx

F1=0;

For j=1:n ,

F1 = F1+2*(a*b*x(i)-f(i))*b*x(i)             

End for

J = [J11, J12, J21, J22]

F = [F1,F2]

Dx=-J/F %Resolver sistema de ecuaciones

X=X+Dx %a=x(1);b=X(2);

XX=[XX,X]

End while

%visualizacion

%Plot()

IV.  Implementación paso a paso de un algoritmo para un modelo no-lineal distinto (variante). Implementado en códigoEditar

____________________________________________________________________________________

I. PRIMERA IMPLEMENTACION:

La primera implementacion esta hecha en base a la aplicacion de un modelo variante

expuesto en la clase del martes 06 de noviembre del 2012.

Este modelo busca continuar experimentando con el trabajo de minimizacion.

El modelo corresponde a un f(t) = A*exp(B*t) con parametros A, B y los datos

(ti,fi), i = 1,2,3, ... n. Se busca formular el algoritmo que minimiza la sumatoria

trabajada anteriormente.

____________________________________________________________________________________

______________________________________________________________________________

SOBRE EL CODIGO:

Este codigo esta explicado en detallo de la misma forma que explicamos al principio

los pasos a seguir para desarrollar un ajuste de datos a un modelo-no lineal.

En el caso anterior fue con el modelo de funcion logistica. Para esta

primera implementacion utilizaremos otro modelo.

El codigo esta implementado en python con una libreria matematica que nos permite

resolver los problemas y que trabaja en terminos muy similares a MATLAB.

SYMPY:

SymPy es una biblioteca escrita en Python cuyo objetivo es reunir todas las caracteristicas de un sistema de algebra computacional (CAS), ser facilmente extensible y mantener el codigo todo lo simple que sea posible. SymPy no requiere ninguna biblioteca externa, salvo para soporte grafico.

Caracteristicas:

En su funcionalidad podemos distinguir entre:

Capacidades basicas, que incluyen:manejo de enteros de precision arbitraria y de numeros racionales,

simplificacion basica, expansion, sustitucion basica,manejo de funciones sobre el cuerpo de los complejos,

derivacion, expansion en series de Taylor o de Laurent,simbolos no conmutativos.

Modulos que incorporan estas tareas:

mas funciones (factorial, zeta, legendre, etc),limites,integracion,divisibilidad y factorizacion de polinomios,

resolucion de ecuaciones algebraicas, diferenciales y sistemas,operaciones con matrices simbolicas,

Algebra de Dirac y de Pauli, Representacion grafica (en 2D y en 3D).

O paquetes externos: symbide: GUI en PyGTK

--> mas informacion sobre sympy en el siguiente link:

http://docs.sympy.org/0.7.2/modules/solvers/solvers.html#

____________________________________________________________________________________

--> CÓDIGO!

#importamos las librerias de sympy... para calculos y variables que

#trabajamos como symbolos

import sympy

from sympy.abc import x, a, b, c

import scipy

from sympy import summation, oo, symbols, log, exp

from sympy import Symbol

from sympy import *

from sympy.plotting import plot

from sympy.plotting import plot3d

#___________________DEFINIMOS EL MODELO (FUNCION f)

A = Symbol('A')

B = Symbol('B')

t = Symbol('t')

f = Function("f")

expresion = A*exp(B*t) #--> MODELO!

f = expresion

#print f

#print  f.diff(t) #para derivar

#print diff(f,t)  #para derivar

#f = A**2

#A=2                                                                                   #probando reemplazos

#f = A**2                                                                          #probando reemplazos

#print "f(2) = " + str(f)   #probando reemplazos

#g = Function("g")

#g = t**2

#print g.subs(t,4).evalf() #para evaluar la funcion f, con respecto a t en un valor cualquiera

#____________#evaluando la funcion f(t), con t = 4

#print f.subs(t,4).evalf()

#__________________________________________________

#_________________DEFINIMOS LA SUMATORIA S

ti = Symbol('ti')

fi = Symbol('fi')

i = Symbol('i') #para iniciar suma

n = Symbol('n') #para termino suma

z = Symbol('z') #funcion a evaluar

#z = (f(ti)-fi)**2 --> expresion de la sumatoria

z  = (f.subs(t,ti).evalf() - fi)**2  #expresion de la sumatoria S, reemplazada en la funcion F

# con expresion = A*exp(B*t) --> f(ti)

#n = 2

#S = summation(z, (i, 1, n)) #sumatoria de S: desde i=1 hasta n

#print S

## PASO1: "Transformar el problema de minimizacion a una tarea de calcular ceros"

# Queremos minimzar S(A,B)

#Calculamos los componentes del vector F:

F1 = Symbol('F1')

F2 = Symbol('F2')

F1 = diff(z,A)

F2 = diff(z,B)

#print F1

#print F2

#Una vez que tenemos F1 y F2 vamos al paso2

## PASO2: "Linealizar el modelo no-lineal"

#Dificultad: Entender y aplicar notacion vectorial

#Entender la secuencia de varios passos de transformacion del problema

#Linealizar el modelo no-lineal, es decir, linealizar

#F(A,B)

#para aproximar vector F(P) = 0, con vector P=(A,B)

#Introducimos Pi= (Ai, Bi), con i=0,1,2...

#formulamos el metodo de newton con la especificacion de F(A,B)

# J11 = Diferencial de F1 con respcto a A

# J12 = Diferencial de F1 con respecto a B

# J21 = Diferencial de F2 con respecto a A

# J22 = Diferencial de F2 con respecto a B

J11 = Symbol('J11')

J12 = Symbol('J12')

J21 = Symbol('J21')

J22 = Symbol('J22')

J11 = diff(F1,A) #term.B

#print J11

J12 = diff(F1,B) #term.AyB

#print J12

J21 = diff(F2,A) #term.AyB

#print J21

J22 = diff(F2,B) #term.AyB

#print J22

ecuaciones_J = [J11,J12,J21,J22]

#print ecuaciones[3]

epsilon = 0.0001

F = [F1, F2]

#resolver sistema de ecuaciones ... usamos from sympy.solvers import solve

equations = [

Eq(-J11, f),

Eq(-J12, f),

Eq(-J21, f),

Eq(-J22, f),

]

#print solve(Eq(-J11, f))

n = 40

valor_ti = 2 #podemso jugar con los parametros

valor_fi = 4 #en distintos valores

result = []

J11 = J11.subs(ti,valor_ti).evalf()

sol1 = summation(J11, (i, 1, n)) #desde i=1 hasta n

aux = solve(-J11)

#result.append( aux )

#print aux

#print J12

J12 = J12.subs(ti,valor_ti).evalf()

J12 = J12.subs(fi,valor_fi).evalf()

sol2 = summation(J12, (i, 1, n)) #desde i=1 hasta n

print sol2

aux = solve(-J12, A) #cambiar si buscamos A o B

result.append( aux[0] )

#print result[0].subs(B,5).evalf() #evaluando valores

J21 = J21.subs(ti,valor_ti).evalf()

J21 = J21.subs(fi,valor_fi).evalf()

sol2 = summation(J21, (i, 1, n))

#print sol2

aux = solve(-J21, A)

result.append( aux[0] )

J22 = J22.subs(ti,valor_ti).evalf()

J22 = J22.subs(fi,valor_fi).evalf()

sol2 = summation(J22, (i, 1, n))

#print sol2

aux = solve(-J22, A)

result.append( aux[0] )

Imagen del código en el editor:

Foto1.jpg

Imagen del código en el editor

            











V. Implementación paso a paso del algoritmo pero en un modelo nuevo (no visto en clases).Editar

Modelo Propio Editar

Si los parámetros del modelo entran en la ecuación en forma no lineal, entonces tenemos un modelo no lineal

¿Lineal o no lineal?

E(Y)=\begin{array}{llll} \beta_0 \ +  \beta_1 \ x_1 +...+ \beta_k \ x_k \qquad\qquad E(Y) = \frac{ \beta_0 \ }{ [1+e^{-( \beta_1 \ + \beta_2 \ x )}]^{ \beta_3 \ }} \end{array}

E(Y)= \begin{array}{llll} \frac{1}{ \beta_0 \ + \beta_1 \ x + \beta_2 \ x^2 } \end{array}

E(Y)= \begin{array}{llll} \beta_0 \ exp(-exp(- \beta_0 \ + \beta_1 \ x)) \qquad E(Y)= \beta_0 \ + \beta_1 \ x+ \beta_2 \ x^2 +...+ \beta_k \ x^k \end{array}

Muchas aplicaciones en biología, medicina,química, agricultura.

En muchos casos surgen a partir de mecanismos físicos, químicos o biológicos conocidos(usando por ejemplo, ecuaciones diferenciales).

Permiten modelar datos reales en forma “natural”con mucha flexibilidad(por ejemplo, con asíntotas,valores positivos de las Y, un único valor máximo).

Muchos parámetros tienen interpretación práctica(tasas de degradabilidad, crecimiento máximo,velocidad máxima de absorción, etc.)

Problemas con los modelos no lineales

Se deben usar métodos iterativos para estimar los parámetros por mínimos cuadrados, máxima verosimilitud, etc.

Es común encontrar problemas de convergencia en los métodos iterativos.

Distintas parametrizaciones pueden afectar no solo la interpretación sino también las propiedades de los estimadores(por ej., una parametrización puede ser muy interesante desde el punto de vista de su interpretación práctica, pero muy mala para lograr convergencia de estimadores).

Familias de funciones no lineales:curvas de crecimiento

Exponencial:  \begin{array}{llll} \frac{dE(Y)}{dt}= \gamma \ E(y);\quad E(Y) = \alpha \ exp( \gamma \ t) \end{array}

Monomolecular:  \begin{array}{llll} \frac{dE(Y)}{dt}= \gamma \ ( \beta \ -Y);\quad E(Y) = \beta_0 \ - \alpha \ exp(- \gamma \ t) \end{array}

Logístico:  \begin{array}{llll} \frac{dE(Y)}{dt}= \gamma \ E(Y)( \beta \ - E(Y));\quad E(Y) = \frac{ \beta \ }{1+exp(- \alpha \ + \gamma \ t)}  \end{array}

Gompertz:  \begin{array}{llll} \frac{dE(Y)}{dt}=  \gamma \ E(Y)[log \beta \ - log E(Y) ];\quad E(Y) = \beta \ exp[- \alpha \ exp(- \gamma \ t ) ] \end{array}

Richards:  \begin{array}{llll} \frac{dE(Y)}{dt}=- \frac{ \gamma \ }{ \kappa \ } E(Y)[( \frac{E(y)}{ \beta \ })^{ \kappa \ } - 1];\quad E(Y) = \beta \ (1 + \kappa \ exp(- \gamma \ (t- \zeta \ )))^{ \frac{-1}{ \kappa \ }} \end{array}


Ejemplos de otras familias de funciones no lineales

Modelo lineal con puntos de cambio:  \begin{array}{llll} E(Y)= \beta_0 \ + \beta_1 \ x I(x \le \alpha \ ) + ( \beta_1 \ a+ \beta_2 \ (x- \alpha \ ))I(x> \alpha \ )    \end{array}


Modelo polinomial recíproco(Holliday):  \begin{array}{llll} E(Y)= ( \alpha \ + \beta \ x + \gamma \ x^2)^{-1} \end{array}

Dentro de los llamados modelos no lineales,  existen 3 tipos;

Reciproco, polinomial y logistico.

Para la Implementacion numero 2, trabajaremos con un modelo Reciproco.


--> CÓDIGO:

Otro modelo NO-Lineal

import sympy

from sympy.abc import x, a, b, c

import scipy

from sympy import summation, oo, symbols, log, exp

from sympy import Symbol

from sympy import *

from sympy.plotting import plot

from sympy.plotting import plot3d

import pylab as P

from matplotlib.transforms import offset_copy

#___________________DEFINIMOS EL MODELO (FUNCION f)

A = Symbol('A')

B = Symbol('B')

Y = Symbol('Y')

x = Symbol('x')

f = Function("f")

expresion = (A+B*x+Y*(x**2))**(-1) #--> MODELO polinomial reciproco (Holliday)

f = expresion

ti = Symbol('ti')

fi = Symbol('fi')

i = Symbol('i') #para iniciar suma

n = Symbol('n') #para termino suma

z = Symbol('z') #funcion a evaluar

z  = (f.subs(Y,ti).evalf() - fi)**2  #expresion de la sumatoria S, reemplazada en la funcion F

# con expresion = A*exp(B*t) --> f(ti)

F1 = Symbol('F1')

F2 = Symbol('F2')

F1 = diff(z,A)

F2 = diff(z,B)

J11 = Symbol('J11')

J12 = Symbol('J12')

J21 = Symbol('J21')

J22 = Symbol('J22')

J11 = diff(F1,A)

J12 = diff(F1,B)

J21 = diff(F2,A)

J22 = diff(F2,B)

ecuaciones_J = [J11,J12,J21,J22]

epsilon = 0.0001

F = [F1, F2]

equations = [

Eq(-J11, f),

Eq(-J12, f),

Eq(-J21, f),

Eq(-J22, f),

]

n = 40

valor_ti = 2 #podemso jugar con los parametros

valor_fi = 4 #en distintos valores

result = []

J11 = J11.subs(ti,valor_ti).evalf()

sol1 = summation(J11, (i, 1, n)) #desde i=1 hasta n

aux = solve(-J11)

J12 = J12.subs(ti,valor_ti).evalf()

J12 = J12.subs(fi,valor_fi).evalf()

sol2 = summation(J12, (i, 1, n)) #desde i=1 hasta n

print sol2

aux = solve(-J12, A) #cambiar si buscamos A o B

result.append( aux[0] )

J21 = J21.subs(ti,valor_ti).evalf()

J21 = J21.subs(fi,valor_fi).evalf()

sol2 = summation(J21, (i, 1, n))

aux = solve(-J21, A)

result.append( aux[0] )

J22 = J22.subs(ti,valor_ti).evalf()

J22 = J22.subs(fi,valor_fi).evalf()

sol2 = summation(J22, (i, 1, n))

aux = solve(-J22, A)

result.append( aux[0] )

#EXPERIMENTOS DE VISUALIZACION

print str(result)

#plot(result[0], (B, -5, 5))

#plot3d(J12, (A, -5, 5), (B, -5, 5)) 

Para la visualizacion vamos a experimentar con la libreria MATPLOTLIB..

esta libreria esta desarrollada para python en base a MATLAB..

los conceptos y sintaxis son muy parecidos...

para mas informacion sobre documentacion y tipos de graficos

en el siguiente link:

http://matplotlib.org/gallery.html

import pylab as P

from matplotlib.transforms import offset_copy

X = P.arange(7)

Y = X**2

fig = P.figure(figsize=(5,5))

ax = P.subplot(2,1,1)

transOffset = offset_copy(ax.transData, fig=fig, x = (result[0].subs(B,7).evalf()), #estamos probando valores

y=(result[0].subs(B,5).evalf()), units='inches')

for x, y in zip(X, Y):

P.plot((x,),(y,), 'ro')

P.show()

Gráfica de resultados:

Foto2.jpg

¡Interferencia de bloqueo de anuncios detectada!


Wikia es un sitio libre de uso que hace dinero de la publicidad. Contamos con una experiencia modificada para los visitantes que utilizan el bloqueo de anuncios

Wikia no es accesible si se han hecho aún más modificaciones. Si se quita el bloqueador de anuncios personalizado, la página cargará como se esperaba.

También en FANDOM

Wiki al azar