Capítulo 7 Aproximación de funciones y ajuste de datos experimentales

22 Pages·2008·0.29 MB·Spanish

Capítulo 7
Aproximación de funciones y ajuste de
datos experimentales
En este capítulo trataremos dos problemas íntimamente ligados. El primero es el problema
delaaproximacióndefuncionesquelopodemosenunciarcomo:
Dada una función f(x) definida en [a,b] y una serie de funciones base ψ (x) definidas tam-
r
biénen[a,b],encontrarloscoeficientesar deformaquelasuma!nr=0arψr(x)sealomáspróxi-
maposiblea f(x)enelintervalo[a,b].
El concepto de proximidad lo definiremos más adelante. El problema de la aproximación es
esencial cuando queremos representar una función en serie de otras más sencillas, como poten-
ciasofuncionestrigonométricas.
El segundo problema surge cuando medimos datos que satisfacen una ley que se comporta
como una función. Típicamente medimos un conjunto de N puntos (x,y), donde la variable
i i
independiente x se supone exacta y todo el error de medida de cada punto se atribuye a la
i
variable dependiente y, que viene afectada de un error experimental σ. Suponemos que la ley
i i
quesatisfacenlosdatossepuededescribirmedianteunmodelodelaformay= f(x)quedepende
de una serie de parámetros a. Nos limitaremos al caso particular en que la dependencia de los
i
parámetroseslineal,esdecir f(x)=!nr=0arψr(x)dondeψr(x)sonfuncionesbaseconvenientes
paradescribirnuestromodeloteóricodelosdatos.Podemosenunciarelsegundoproblemacomo:
Determinarlosvaloresdelosparámetrosa quehacenquelacantidad
i
χ2(a ,a ,...a )= !N (yi−!nr=0arψr(xi))2
0 1 n σ2
i=1 i
seamínima.
Esteeselproblemadelmodeladodedatosexperimentales.Ambosproblemas,aproximación
de funciones y modelado de datos, están íntimamente ligados y comparten las mismas técnicas
deresolución.
115
116CAPÍTULO7. APROXIMACIÓNDEFUNCIONESYAJUSTEDEDATOSEXPERIMENTALES
7.1. Proximidad de funciones: Distancias y Normas
Enprimerlugar,hayquedefinirelconceptodeproximidaddedosfuncionesenunintervalo.
Paraellohayqueintroducirunadistanciaentrelasdosfunciones.Lasdistanciassesuelendefinir
mediante normas. Si tenemos una norma definida para funciones f(x) , se define la distancia
" "
entre dos funciones f(x) y g(x) como d(f(x),g(x))= f(x) g(x) . Hay diversas normas utili-
" − "
zadasfrecuentemente.LamásutilizadaeslanormademínimoscuadradosoL definidacomo
2
b
2
f(x) g(x) = (f(x) g(x)) dx
" − "2 −
a
!
enunintervaloycomo
n
f(x) g(x) = !(f(x) g(x))2
" − "2 i − i
i=0
sobreunconjuntodiscretodepuntos.EngenerallanormaL sedefinecomo
p
b
p
f(x) g(x) = f(x) g(x) dx
" − "p | − |
a
!
sobreunintervaloycomo
N
f(x) g(x) = ! f(x) g(x) p
" − "p | i − i |
i=1
sobre un conjunto discreto de puntos. En aproximación de funciones, aparte de la norma L , se
2
utilizanusualmentelanormaL ylallamadanormaL ,definidacomo
1 "
f(x) g(x) =ma´x f(x) g(x)
" − "" | − |
sobre un intervalo o conjunto discreto de puntos. La aproximación de funciones que minimiza
la norma L se conoce como aproximación minimax. Cuando deseamos una aproximación a
"
una función en un intervalo por otra más sencilla, la aproximación minimax es quizás la más
razonable, ya que limita el error máximo cometido en un punto arbitrario del intervalo. Sin
embargo, cuando tenemos puntos experimentales afectados de un error estadístico, entonces la
aproximación de mínimos cuadrados, en la versión de mínimo χ2, es la única justificada desde
elpuntodevistaestadístico.
7.2. Aproximación de mínimos cuadrados
7.2.1. Normas a partir de productos escalares
Sidefinimoselproductoescalardedosfuncionescomo
b
< f(x) g(x)>= f(x)g(x)dx
|
a
!
7.2. APROXIMACIÓNDEMÍNIMOSCUADRADOS 117
sobreunintervaloy
N
< f(x) g(x)>= ! f(x)g(x)
i i
|
i=1
sobre un conjunto discreto de puntos. La norma L se puede escribir en función del producto
2
escalarcomo
f(x) g(x) =< f(x) g(x) f(x) g(x)>
" − "2 − | −
tantosobreunintervalocomounconjuntodiscretodepuntos.
7.2.2. Las ecuaciones normales de mínimos cuadrados
Engeneraldeseamosaproximarunafunción f(x)porunacombinaciónlinealdeunconjunto
den+1funcionesbaseψ (x)
r
n
f(x)= !a ψ (x)
r r
r=0
El caso más frecuente es cuando ψ (x)=xr , que se denomina aproximación polinómica. Para
r
lleva a cabo la aproximación tenemos que encontrar los coeficientes a , a ,...,a que hacen la
0 1 n
función
n
E(a ,a ,...,a )= f(x) !a ψ (x)
0 1 n r r
" − "
r=0
" "
" "
mínimo.TenemosqueminimizarE considerada"comounafunción"delosparámetrosa ,
r
" "
n n
E(a ,a ,...,a ) = < f(x) !a ψ (x) f(x) !a ψ (x)>=
0 1 n r r r r
− | −
r=0 r=0
n n
< f(x) f(x)> 2!a < f(x) ψ (x)>+ ! a a <ψ (x) ψ (x)>
r r r s s r
| − | |
r=0 r,s=0
Las condiciones que se deben de cumplir para que exista un mínimo son, en primer lugar, la
anulación de las derivadas primeras con respecto de los parámetros, y en segundo lugar que la
matrizdederivadassegundasoHessianoseadefinidapositiva
∂E(a ,a ,...,a )
0 1 n
= 0
∂a
i
∂2E(a ,a ,...,a )
0 1 n
> 0
∂a∂a
# i j #
# #
# #
Laprimeradelascondiciones#da #
∂E(a ,a ,...,a ) n
0 1 n
= 2< f(x) ψ(x)>+2!a <ψ (x) ψ(x)>=0
i r r i
∂a − | |
i r=0
118CAPÍTULO7. APROXIMACIÓNDEFUNCIONESYAJUSTEDEDATOSEXPERIMENTALES
Estacondiciónimplicaelcumplimientodeunsistemadeecuaciones
n
!a <ψ (x) ψ(x)>=< f(x) ψ(x)> (7.1)
r r i i
| |
r=0
queseconocencomoecuacionesnormales.Constituyenunsistemalinealparalosparámetros
Aa=b
dondeaeselvectordeparámetros,belvectordetérminosindependientesyAlamatrizdecoefi-
cientes.Lasegundacondiciónsecumplesiempre,loquesepuedeverexplícitamentesuponiendo
quevariamoslosparámetrosa a +δa ycalculamosladiferencia
r r r
→
E(a +δa ,a +δa ,...,a +δa ) E(a ,a ,...,a )=
0 0 1 1 n n 0 1 n
−
n n
< f(x) !(a +δa )ψ (x) f(x) !(a +δa )ψ (x)>
r r r r r r
− | −
r=0 r=0
n n
< f(x) !a ψ (x) f(x) !a ψ (x)> =
r r r r
− − | −
r=0 r=0
n n n n
= 2!δa <ψ (x) f(x) !a ψ (x)>+< !δa ψ (x) !δa ψ (x)>
r r s s r r r s
− | − |
r=0 s=0 r=0 r=0
El primer término se anula por el cumplimiento de las ecuaciones normales y el segundo es
estrictamentepositivo,puestoqueeslanormadeunvectornonulo.
Elcasomássimpleescuandotenemosúnicamentedosfuncionesbaseψ yψ .Entonceslas
0 1
ecuacionesnormalesquedancomo
a < ψ ψ >+a <ψ ψ >=<ψ f >
0 0 0 1 0 1 0
| | |
a < ψ ψ >+a <ψ ψ >=<ψ f >
0 1 0 1 1 1 1
| | |
cuyassoluciones,aplicandolafórmuladeCramerson
<ψ f > <ψ ψ >
0 1 0
| |
<ψ f > <ψ ψ >
1 1 1
a = # | | #
0 #<ψ ψ > <ψ ψ >#
# 0 0 0 1 #
| |
#<ψ ψ > <ψ ψ >#
1 0 1 1
# | | #
# #
<ψ ψ > <ψ f >
# 0 0 0 #
# | | #
<ψ ψ > <ψ f >
1 0 1
a = # | | #
1 #<ψ ψ > <ψ ψ >#
# 0 0 0 1 #
| |
#<ψ ψ > <ψ ψ >#
1 0 1 1
# | | #
# #
# #
# #
7.2. APROXIMACIÓNDEMÍNIMOSCUADRADOS 119
Siconsideramoselcasodelajustelineal,ψ =1yψ =x,enelcasodeunconjuntodiscretode
0 1
puntostenemos
N N
<ψ0|ψ0 >= !1=N, <ψ0|ψ1 >=!Ni=1xi, <ψ1|ψ1 >= !xi2,
i=1 i=1
N
<ψ0|f >= ! f(xi) <ψ1|f >=!Ni=1xif(xi)
i=1
Poniendo y = f(x) tenemos las fórmulas usuales del ajuste de un conjunto de puntos por míni-
i i
moscuadrados:
a = !Ni=1yi!Ni=1xi2−!Ni=1xi!Ni=1xiyi a = !Ni=1yi!Ni=1xi2−N!Ni=1xiyi
0 N!Ni=1xi2− !Ni=1xi 2 1 N!Ni=1xi2− !Ni=1xi 2
$ % $ %
En el caso de aproximaciones polinómicas de orden más elevado (parabólicas, cúbicas, o com-
binaciones lineales de varias potencias distintas) procederíamos de forma análoga, resolviendo
las ecuaciones por uno de los métodos vistos en el capítulo 4, en vez de por la regla de Cramer.
Podemos pensar que podemos continuar de esta forma hasta cualquier orden de aproximación
aunque este no es el caso. De hecho para más de 10 funciones, las ecuaciones normales están
mal condicionadas, y dan resultados imprecisos con doble precisión. Para orden 100, incluso
con cuádruple precisión en procesadores de 64 bits se obtienen resultados muy imprecisos. Sin
embargo no es raro que sea necesario aproximar una función por varios centenares de funciones
base.Estoocurreporejemplocuandosedescomponeunaondasonoraenarmónicosocuandose
estudian imágenes. Si obtenemos una solución imprecisa de las ecuaciones normales los agudos
de una onda serían incorrectos y la imagen no sería nítida. Por ello hace falta un método eficaz
deevitarelmalcondicionamiento.Elloseconsigueconfuncionesortogonales.Decimosquelas
funcionesψ sonortogonalessi
r
<ψ ψ >=n δ
r s r rs
|
donde n es la normalización de la función y δ es la delta de Kronecker. En este caso las
r ij
ecuacionesnormalessesimplificana
a <ψ ψ >=<ψ f >
r r r r
| |
conlasolución
<ψ f >
r
a = |
r
<ψ ψ >
r r
|
La utilización de funciones ortogonales tiene dos ventajas: la primera es que desaparece el mal
condicionamiento, y la segunda es que cada coeficiente es independiente de los demás. Por lo
tanto, si deseamos extender la aproximación a un orden superior, los coeficientes ya calculados
no varían, por lo se dice que tienen la propiedad de permanencia. Esta independencia es muy
importante en el caso de datos experimentales, puesto que implica que los distitos coeficientes
obtenidosajustandomediantefuncionesortogonalesnoestancorrelacionadosesdadísticamente.
120CAPÍTULO7. APROXIMACIÓNDEFUNCIONESYAJUSTEDEDATOSEXPERIMENTALES
7.2.3. Series de Fourier
Sin duda alguna, las funciones ortogonales más utilizadas son las funciones trigonométricas
sin(x) y cos(x). El conjunto de funciones 1,cos(x),sin(x),cos(2x),... son ortogonales en el
{ }
intervalo[ π,π]conlasrelacionesdeortogonalidad
−
π π π
dxcoskxcosmx= dxcoskxsinmx= dxsinkxsinmx=0 m=k
$
π π π
!− !− !−
π π
dxcoskx= dxsinkx=0 k>0
π π
!− !−
π π π
dx(coskx)2 = dx(sinkx)2 =π dx=2π
π π π
!− !− !−
Eldesarrollodeunafuncióncomo
a "
0
f(x) + !(a cosrx+b sinrx)
r r
∼ 2
r=1
se conoce como serie de Fourier. Converge en la norma de mínimos cuadrados siempre que la
función sea periódica en [ π,π] y continua. Cuando la serie se trunca a un número finito de
−
términos, frecuentemente grande, tenemos la aproximación de Fourier. Los coeficientes vienen
dadospor
1 π 1 π 1 π
a = dxf(x) a = dxf(x)cosrx b = dxf(x)sinrx
0 r r
π π π
π π π
!− !− !−
En casos analíticamente sencillos los coeficientes de Fourier se calculan fácilmente. Consi-
deremosporejemplounaondacuadrada,queseutilizafrecuentementeenelectrónica.
1 π x<0
f(x)= − − ≤
1 0 x<π
& ≤
Esta función es una función impar. También es discontinua, pero a pesar de esto la serie de
Fourier converge. Como cosx es par, los coeficientes a se anulan. Los coeficientes b vienen
r r
dadospor
1 π 2 π 2 π 0 rpar
b = dxf(x)sinrx= dxsinrx= cosrx = 4
r
π !−π π !0 π #0 ’ πr rimpar
#
2 " sin[(2r+1)x] #
#
f(x) !
∼ π 2r+1
r=0
EnelcasodeunafunciónperiódicadeperíodoT,eldesarrollotomalaforma
a " 2πrt 2πrt
0
f(x) + !(a cos +b sin )
r r
∼ 2 T T
r=1
con
2 T/2 2 T/2 2πrt 2 T/2 2πrt
a = dtf(t) a = dtf(t)cos b = dtf(t)sin (7.2)
0 r r
T T T T T
T/2 T/2 T/2
!− !− !−
7.2. APROXIMACIÓNDEMÍNIMOSCUADRADOS 121
SeriedeFourierdiscreta
Las funciones trigonométricas también son ortogonales sobre un conjunto finito de puntos.
Dada una función f(t) periódica con período T, si tomamos un conjunto de N+1 puntos igual-
menteespaciadosentre0yT (t =sT/(N+1),s=0,...,N)sesatisfacenlassiguientesrelacio-
s
nesdeortogonalidad
2πkt 2πmt N 2πks 2πms 0 k=m,k=m=0,N+1
<sin sin >= !sin sin = $
N+1
T | T N+1 N+1 k=m=0,N+1
s=0 ’ 2 $
2πkt 2πmt N 2πks 2πms
<sin cos >= !sin cos = 0
T | T N+1 N+1
s=0
0 k=m
2πkt 2πmt N 2πks 2πms $
N+1
<cos cos >= !cos cos = k=m=0,N+1
T | T N+1 N+1  2 $
s=0  N+1 k=m=0,N+1

Eldesarrollo
a n 2πks 2πks
0
f(t) + ! a cos +b sin
k k
∼ 2 N+1 N+1
k=1+ ,
converge a f(t) sobre el conjunto de N+1 puntos en el sentido de mínimos cuadrados. Cuanto
tomamos N+1 coeficientes, el desarrollo interpola a la función f(t) en el conjunto de N+1
puntos.SiN espar(númerodepuntosimpar),lafuncióninterpoladoraes
N/2
sT a 2πks 2πks
0
F = + ! a cos +b sin
N+1 k k
N+1 2 N+1 N+1
+ , k=1+ ,
mientrasquesiN esimpar(númeropardepuntos)
sT a0 (N−1)/2 2πks 2πks a(N+1)/2
F ( )= + ! a cos +b sin + cosπs
N+1 k k
N+1 2 N+1 N+1 2
k=1 + ,
Loscoeficientesdeldesarrollovienendadospor
2 N sT 2πks
a = ! f cos
k
N+1 N+1 N+1
s=0 + ,
2 N sT 2πks
b = ! f sin
k
N+1 N+1 N+1
s=0 + ,
Es interesante notar que a y b vienen dados por la evaluación numérica mediante la regla
k k
trapezoidal para N+1 intervalos (N+2 puntos, ampliando con el extremo del t = T) de las
integrales de las ecuaciones 7.2, notando que f(0) = f(T), T = (N+1)h, y que los senos se
122CAPÍTULO7. APROXIMACIÓNDEFUNCIONESYAJUSTEDEDATOSEXPERIMENTALES
anulanenlosextremosdelintervalo:
2 f(0)+ f(T) N sT 2πks
a = +! f cos
k
N+1 2 N+1 N+1
- s=1 + , .
2 N sT 2πks
b = ! f sin
k
N+1 N+1 N+1
-s=1 + , .
7.3. Polinomios ortogonales
El conjunto más sencillo de funciones ortogonales son los polinomios. Se pueden definir
sobreunconjuntodiscretodepuntososobreunintervalocontinuo.Vamosadefinirlosporahora
con coeficiente de la potencia más elevada igual a la unidad. De esta forma siempre existe una
relaciónderecurrenciadeltipo(p (x) xp (x))
k+1 k
−
k
p (x)=xp (x)+!ck+1p (x) (7.3)
k+1 k s s
s=0
ya que (p (x) xp (x)) es un polinomio de grado k, y por lo tanto siempre se puede expresar
k+1 k
−
como combinación lineal de p (x),...,p (x). Vamos a suponer únicamente la existencia de un
0 k
producto escalar sobre un intervalo [a,b] o sobre un conjunto discreto de N+1 puntos. Dicho
productoescalarlossupondremosdelaformamásgeneralconunafunciónpesow(x)enelcaso
continuoyunconjuntodepesosw enelcasodiscreto
s
b
dxw(x)p (x)p (x)
< p (x) p (x)>= a k j
k | j & /!Ns=1wspk(xs)pj(xs) 0
Tenemos que determinar los coeficientes ck+1. Para ello multiplicamos escalarmente la ec. 7.3
s
porunpolinomiodado p (x),r k,
r
≤
k
< p p >=0=< p xp >+!ck+1 < p p >=< p xp >+ck+1 < p p >
r k+1 r k s r s r k r r r
| | | | |
s=0
dedonde
< p xp >
ck+1 = r| k
r
− < p p >
r r
|
Como < p xp >=< p x p > y xp (x) es un polinomio de grado r+1, que se puede expresar
r k r k r
| |
comounacombinaciónlinealde p ,...,p ,< p xp >=0parar=0,1,...,k 2.Porlotanto
0 r+1 r k
| −
sólock+1 yck+1 puedenserdistintosde0.Vienendadospor
k 1 k
−
< p xp >
ck+1 = k−1| k
k 1 −< p p >
− k 1 k 1
− | −
7.3. POLINOMIOSORTOGONALES 123
y
< p xp >
ck+1 = k| k
k − < p p >
k k
|
Lospolinomiosortogonalessatisfacenporlotantolarelaciónderecurrencia
p (x)=(x+ck+1)p (x)+ck+1p (x)
k+1 k k k 1 k 1
− −
Paraqueestarelaciónsecumplatambiénpara p (x)sedefine p (x)=0.Paraobtenerelajuste
1 1
−
por mínimos cuadrados de una función dada f(x) , sólo tenemos que calcular los coeficientes
ck+1 yck+1 mediantelasecuacionesanterioresparaobtenerlospolinomiosnecesariosmediante
k k 1
larelació−nderecurrencia.Elajustedemínimoscuadradosdeordennvienedadopor
n
!a p (x)
r r
r=0
dondea seobtienede
r
< f p >
r
a = |
r
< p p >
r r
|
El incremento del orden de aproximación en una unidad implica, por lo tanto, el cálculo de un
nuevopolinomioyuncoeficiente,loqueequivalearealizar6productosescalares,quesereducen
a4dadolasconstantesdenormalizacióndelospolinomios< p p >sehancalculadodurantela
r r
|
obtencióndelcoeficienteprevio.Estaeslaformamáseficientedeajustardatosmediantepolino-
miosdeordenelevado,tantoparadatosdiscretoscomocontinuos,puesseevitanerroresdebidos
al mal condicionamiento de las ecuaciones normales, y por otro lado el esfuerzo numérico es
menor, y se puede elevar el orden aprovechando los cálculos realizados para un orden inferior.
En el caso de datos discretos, el único inconveniente es la dependencia de los polinomios del
conjunto de puntos, lo cual no es importante, pues la suma de polinomios ortogonales se puede
expresardeformainmediatacomounpolinomioordinario.
Paradatosdefinidosenintervaloscontinuoshaypolinomiosortogonalesbienconocidospara
diversospesoseintervalos,algunosdeloscualessedanenlatabla7.1
Tabla7.1:Principalespolinomiosortogonales
Nombre Peso Intervalo Símbolo
Legendre 1 [ 1,1] P (x)
n
−
Hermite exp( x) [ ","] H (x)
n
− −
Laguerre exp( x2) [0,"] L (x)
n
−
Chebychev 1/√1 x2 [ 1,1] T (x)
n
− −
Chebychev2ªespecie √1 x2 [ 1,1] U (x)
n
− −
Silafunción f seconoceanalíticamenteosepuedecalcularconfacilidadencualquierpunto
que se desee, los coeficientes del desarrollo de la función en serie de polinomios ortogonales se
puedencalcularporcualquieradelosmétodosdeintegraciónvistosenelcapítuloanterior.
124CAPÍTULO7. APROXIMACIÓNDEFUNCIONESYAJUSTEDEDATOSEXPERIMENTALES
7.3.1. Serie de Chebychev discreta
Otro conjunto de funciones que satisfacen relaciones de ortogonalidad sobre un conjunto
discretodepuntossonlospolinomiosdeChebychev.
7.4. Aproximación minimax
7.5. Aproximación por funciones racionales
7.6. Modelado de datos experimentales
7.6.1. Variables aleatorias, valores esperados y varianzas
Una variable aleatoria es una variable que puede tomar un conjunto de valores (continuo o
discreto) y que cada valor aparece con una probabilidad determinada. Por ejemplo el valor de la
cara de un dado puede tomar 6 valores con probabilidad 1/6. El número de desintegraciones de
una muestra radioactiva en la unidad de tiempo toma valores enteros. La variable puede tomar
valorescontinuos,encuyo casoexisteunadistribucióndeprobabilidado densidaddeprobabili-
dad p(x), definida en [ ",+"]. La probabilidad de que x tome un valor comprendido entre dos
−
valoresaybvienedadapor
b
P(a<x<b)= p(x)dx
a
!
Sedefineelvaloresperadodex,E[x],tambiéndenominadovalormedio,como
"
E[x]=x= xp(x)dx
"
!−
ylavarianzaσ2(x)como
"
σ2(x)=E[(x x)2]= (x x)2p(x)dx
− −
"
!−
Frecuentemente tenemos varias variables aleatorias que pueden aparecer simultáneamente.
En este caso tenemos una distribución de probabilidad conjunta p(x ,x ,...,x ). Si tenemos dos
1 2 n
variablesaleatoriasx yx ,sedefinelacovarianzaσ(x ,x )como
1 2 1 2
"
σ(x ,x )=E[(x x )(x x ]= (x x )(x x¯ )p(x ,x )dx dx
1 2 1 1 2 2 1 1 2 2 1 2 1 2
− − − −
"
!−
Si dos variables son independientes, su covarianza se anula, ya que en este caso p(x ,x ) =
1 2
p(x )p(x ) y la integral anterior se descompone en el producto de dos integrales que se anulan,
1 2
locualsedemuestrafácilmenteteniendoencuentaladefinicióndelvalormedio.
Los datos experimentales se comportan como variables aleatorias. Cada vez que medimos
unamagnitudfísicaconsuficienteprecisiónobtenemosunvalordistinto.Elconjuntodevalores
deunaseriedemedidassedistribuyeconunafuncióndedistribucióndeprobabilidad.Unaserie
demedidasx secaracterizaporsuvalormediox¯ysudesviacióntípicaσ .
i x