Chapra17-18

5/14/2018 Chapra17-18 - slidepdf.com

http://slidepdf.com/reader/full/chapra17-18 1/74

C A P I T U L O 1 7Regres ión por mínimoscuadrados

Donde se asocian errores sustanciales con los datos, la interpolación polinomial esinapropiada y pue de dar resultados insatisfactorios cuando se usa para predecir valor es

intermedios. Por ejemplo, en la f igura 17.1a se muestran siete datos derivados experi-mentalmente que exhiben variabilidad significativa. Una inspección visual de dichosdatos sugiere una posible relación entre y y x. Es decir, la tendencia general indica quelos valores más altos de y son asociados con los valores más altos de x. Ahora, si unainterpolación de sexto orden se ajusta a estos datos (figura 17.1¿»), pasará justo a travésde todos los p untos. Sin embargo, a causa de la variabilidad en los datos, la curva o scilaen forma amplia en el intervalo entre los puntos. En particular, los valores interpoladosen x = 1.5 y JC = 6.5 parecen estar mu y adelante del rango suger ido por los datos.

Una estrategia más apropiada para tales casos es derivar una función aproximadaque ajuste la forma d e la tendencia general de lo s datos sin ajustar necesariamente conlos puntos individuales. La figura 17.1c ilustra cóm o se pued e usar por lo general unalínea recta para caracterizar la tendencia de los datos sin pasar a través de un punto en

particular.Una manera para determinar la línea en la f igura 17.1c es inspeccionar en forma

visual los datos graneados y después trazar una "mejor" línea a través de los puntos.Aunque tales procedimientos por "vistazo" apelan al sentido común y son válidos paracálcu los superficiales, resultan deficientes por ser arbitrarios. Es decir, a meno s que lospuntos definan una línea recta perfecta (en tal caso la interpolación podría ser apropiada), diferentes analistas podrían dibujar distintas líneas.

Para hacer a un lado la subjetividad se debe concebir algunos cr it er is j^ n el f inde establecer una base para el ajuste. Una forma de hacerlo es derivar una curva queminimice la discrepancia entre los puntos y la curva. Una técnica para cumplir con talobjet ivo se conoce como regresión por mínimos cuadrados, que se analizará en estecapítulo.

1 7 . 1 R E G R E S I Ó N L I N E A L

El ejemplo más simple de una aproximación por mínimos cuadrados es mediante elajuste de un conjunto de pares de observ aciones: (x l,y l), (x 2,y 2), • • •. (x„ ,y„) a una línearecta. La ex presión m atemática para esta última es

V = c/<) + + (' (17 .1)



466 REGRESIÓN POR MÍNIMOS CUADRADOS

a)

b)

F I G U R A 1 7 . 1a) Datos que exhiben unerror significativo, b) Ajustepolinomial oscilando másallá del rango de datos,c) Resultados mássatisfactorios mediante elajuste por mínimoscuadrados. c)

donde a 0 y a¡ son coef icientes que representan el in tercepto y la pendiente, respectivamen te , y e es el er ror, o res iduo, entre el modelo y las observacione s, las cuales se pueden representar al reordenar la ecuación (17 .1) como

e = v - ero - a\x

Así, el envr o residuo es la discrepancia entre el valor real <icy y el valor aproximado, o „

+ predicho por la ecuación lineal.



1 7 ,1 R E G R E S I Ó N L I N E A L 4 * 9

1 7 . 1 . 1 C r i t e ri o s p a r a u n " m e j o r " a j u st eUna estrategia para ajustar a la ¿"mejor"? l ínea a través de los datos podr ía ser minim izarla suma de los er rores res iduales para todos los datos d isponibles , como en

J 2e

' <= (-' " °° " a ix

> ) (17.2)

í = l 1 = 1

donde n — núm ero to tal de puntos . Sin embargo, és te es un cr i ter io inadecua do, com o loilustra la f igura 17 .2a, la cual muestra el ajuste de una l ínea rec ta de dos pun tos . O bvia-

F I G U R A 1 7 . 2Ejemplos de algunos criterios para "el mejor ajuste" inadecuados para regresión:a) minimiza la suma de los residuos, b) minimiza la suma de los valores absolutos de losresiduos y c) minimiza el error máximo de cualquier punto individual

xc)



4 68 R E G R E S I Ó N POR M Í N I M O S C U A D R A D O S

m e n t e , el mejor ajuste es la línea que conec ta los pun tos . Sin embargo, cualquier l ínearecta que pasa a t ravés del pun to med io que conec ta la l ínea (excepto una l ínea perfectaver t ical) resulta en un va lo r mín imo de la ecuación (17 .2) igual a cero debido a los

errores que se cance lan .Por tanto , o tro cr i ter io lógico podr ía ser min imizar la s u m a de los valores absolu tos

de las d iscr epanc ias , com o en

n n

^ k l = Yl ^ y> - a o - a i x í \1=1 í = i

La figura 17.2b demues t r a por qué este cr i ter io es t amb ién inadecuado . Para los cua t ropuntos ex puestos , cualquier l ínea recta que esté dentro de l a s l íneas pun teadas min imiza

rá el valor absolu to de la suma. Así , es te cr i terio tampo co da un único mejor ajuste.Una tercera es trategia para ajustar a la mejor l ínea es el cr i ter io mirúmax. E

técnica, la l ínea se elige de manera que m i n i m i c e la máxima d is tanc ia que t enga un

punto indiv idual desde la l ínea. Como se i lustra en la figura 17.2c, tal estrategia no es

adecuada para regresión , ya que t iene una excesiva inf luencia en puntos fuera del conj u n t o ; es decir, un solo punto con un gran er ror . Deber ía observarse que el pr incip iom i n i m a x es algunas ocasiones muy adecuado para ajustar una simple función a una

complicada función (Carnahan, Luther y Wilkes, 1969) .Una estrategia que supera los defectos de los proced imien tos mencionados es min i

m i z a r la suma de los cuad rados de los residuos entre la y m e d i d a y la y ca lcu lada con el

modelo l ineal

n n n^ = 2 ^ = 2 ( y ; , m e d i d a - > ' 1 , m o d e l o ) 2= 2<>,-ao -«l^) 2

( 1 7 - 3 )

k ;=1 i=l ¡=1

Este cr i ter io tiene va r ias ventajas, en tre el las el hecho de que se ob t iene una l ínea únicapara un cier to conjunto de datos . Ante s de analizar esas propiedades, presentaremos una

técnica para determinar los valores de a0 y a{ que m i n i m i z a n la ecuac ión (17 .3 ) .1 7 . 1 . 2 A j u s t e por m í n i m o s c u a d r a d o s de una l í ne a re c t a

Para de te rminar los valores de a 0 y av la ecuación (17 .3) es diferenciada con respecto a

cada coef iciente:

2Yl(y¡ ~a° ~ aiXi)

2

^ 2 [ { y¡

- « o - « i )• * • ; ]

Observe que hem os s imp l i f icado los s ímbo los de la sumator ia; a m e n o s que se indiqueotra cosa, todas las sumator ias son de i = 1 hasta n. Al f i jar esas de r ivadas igual a cero ,resultará en un m í n i m o S

r. Si se hace esto , las ecuaciones se pueden exp resar como

3S r

da0

dS r9<2l



17.1 REGRESIÓN LINEAL 4 1 9

Ahora, si hacemos que £a 0 = na0, podemos expresar las ecuaciones como un conjuntode dos ecuaciones lineales con dos incógnitas (a 0 y a,):

na0+ (X!') i ~ X̂ ' O 7 ' 4 )a 0 + Q T X ? ) a i = £ > y ¡ . (17.5)

Estas son llamadas ecuaciones normales, y pueden ser resueltas en forma simultánea

(17.6)nlx* - (Zx,.)*

Este resultado, entonces, se puede usar en conjunto con la ecuación (17.4) y resolverpara

a0 = y - axx (17.7)

donde y y x son las medias de y y x, respectivamente.

E J E M P L O 17 .1 Regresión lineal

Enunciado del problema. Ajuste a una línea recta los valores de x y y en las dos

i primeras columnas de la tabla 17.1.

Solución. Se calculan las siguientes cantidades:

n = ix

>-y¡ =1 1 9 5

J2x

<2 = 1 4 0

^ J C , - = : 2 8 i = y = 4

^ 24y¡ = 24 x = — = 3.428571

Mediante las ecuaciones (17.6) y (17.7)

7(119 .5) -28(24) „„„„„„„ai = — -V- = 0.8392857

7(140) - (28) 2

a 0 = 3.428571 - 0.8392857(4) = 0.07142857

T A B L A 17.1 Cálculos para un análisis de error del ajuste l ineal.

_*f y> (y;-y)2 (yz-op-*»!*;)2

1 0.5 8.5765 0.16872 2.5 0.8622 0.5625

1 3 2.0 2.0408 0.3473j 4 4.0 0.3265 0.3265í 5 3.5 0.0051 0.5896I 6 6.0 6.6122 0.7972

7 ¿5 . 54,29Qg 01 99 3S 24.0 22.7143 2.9911



4 7 0 R E G R E S I Ó N P O R M Í N I M O S C U A D R A D O S

Por tanto , el ajuste por mínimos cuadrados es

J y = 0 .07142857 + 0 .8392857*

S La l ínea, jun to con los datos , son mos trados en la figura 17 .1c.

1 7 . 1 . 3 C u a n t if ic a c i ó n d e l e r r o r d e u n a r e g r e s i ó n l i n e a l

Cualquier o tra l ínea que la calculada en el ejemplo 17.1 resulta en una gran suma decuadrados de los res iduos. Así , la l ínea es única y en término s de nuestro cr i ter io elegidoes una l ínea "mej or" a través de los punt os . Un número adicional de propied ades de e steajuste se puede elucidar al examinar más de cerca la forma en que se calcularon losresiduos. Recuerd e que la suma de los cuadrad os se def ine com o [véase ecuació n (17 .3)]

Observe la s imili tud entre las ecuaciones (PT5.3) y (17 .8) . En el pr imer caso , elcuadrado del res idu o representa el cuadrado d e la d iscrepancia entre los datos y una solaest imación de la medida de tendencia central ( la media) . En la ecuación (17 .8) , el cuadrado de los res iduos representa el cuadrado de la d is tancia ver t ical entre los datos y o tramed ida de tendenc ia central : la l ínea recta (véase f igura 17 .3) .

La analogía se puede extender más para casos donde 1) la d ispers ión de los puntos

alrededor de la l ínea es de magnitud s imilar junto con todo el rango de datos , y 2) ladis tr ibución de esos puntos cerca de la línea es normal. Se pued e dem ostrar que s i es toscr i ter ios se cumplen, la regresión por mínimos cuadrados proporcionará la mejor (esdecir , una de las mejores) es t imación d e a 0 y a { (Draper y Smith , 1981) . Esto es cono cidoen estadís t ica como el principio de probabilidad máxim a. Además, s i es tos cr i ter ios se

El residuo en la regresión lineal representa la distancia vertical entre un dato y la línea recta.

(17 .8)í = i

F I G U R A 1 7 . 3

Medición

a 0 + a,*,

x



17 . 1 REGRE S I ÓN L I NEAL

cumplen , una "desviac ión estándar" para la l ínea de regresión se puede dete rmina r c om o[compare con la ecuación (PT5.2)]

(17 .9)

donde sy/x es l lamado el error estándar del estimado. La no tac ión de l sub índ ice "y/x "

designa que el er ror es para un valor predich o de y correspondiente a un valor par t icularde x. También, observe que ahora d iv id imos entre n — 2 debido a los dos datos es t imados (a Q y Í Z J ) , que se usaron para calcular S r; así , se t iene dos grados d e l iber tad . Com o lohicimos en nuestro anális is para la desviación estándar en P T 5 . 2 . 1 , otra just if icaciónpara dividir entre n — 2 es que no exis te algo como "d atos d isp ersos " alrededor de u nalínea recta que conecte dos puntos . De esta manera, para el caso donde n — 2, la ecuación (17 .9) da un resultado s in sentido al inf in ito .

Justo como fue el caso con la desviación estándar , el er ror es tándar de la es t ima cióncuantif ica la d ispers ión de los datos . Sin embargo, s y/x cuantif ica la d ispers ión alrededor

de la línea de regresión, com o se mue stra en la figura 17 .4¿, en contraste con la desviación estándar or ig inal s y que cuantifica la dispersión alrededor de la media ( f igura 17 .4a) .

Los conceptos anter iores se pueden usar para cuantif icar la "bondad" de nuestroajuste. Esto es en par t icular ú t i l para compa rar d iferentes regresiones (véase f igura 17 .5) .Para hacer es to , regresamos a los datos or ig inales y determinamos la suma total de los

cuadrados alrededor de la me dia para la var iable dependiente (en nuestro caso , y). C o m ofue el caso para la ecua ción (PT 5.3) , es ta cantidad se designa por S¡ . Ésta es la magnituddel er ror res idual asociado c on la var iable dependiente antes de la regresión . Despué s derealizar la regresión , calculamos 5 r , la suma de los cuadrados de los res iduos alrededor

de la l ínea de regresión . Esto caracter iza el er ror res idual que qu eda de spué s de la regres ión . Esto es , por tanto , algunas veces l lama do la suma inexplicable d e los cuadrados. La

F I G U R A 1 7 . 4

Datos de regresión que muestran o) la dispersión de los datos alrededor de la media de la variable dependiente y b) ladispersión de los datos alrededor de la mejor línea de ajuste. La reducción en la dispersión va de a) a b) , como lo indicanlas curvas en forma de campana a la derecha, representan la mejora debida a la regresión lineal.




F I G U R A 1 7 . 5Ejemplos de regresión lineal con errores residuales a) pequeños y b) grandes.

diferencia entre estas dos cantidades, S t — S r, cuantifica la mejora o reducción de errordebido a que describe lo s datos en términos de una línea recta en vez de c om o un valorprome dio. Co mo la magnitud de esta cantidad es dependiente de la escala, la diferenciaes normalizada a S t para obtener

r2

=S t ~ S r ( 1 7 . 1 0 )

S,donde r 2 es conocido como el coeficiente de determinación y r es el coeficiente de conf

lación (— Vr5). Para un ajuste perfec to, S r — 0yr = r 2 = 1, significa que la líneaexplica el 100% de la variabilidad de los datos. Para r = r 2 = 0 , S r = S, y el ajuste norepresenta ninguna mejora. Una formulación alternativa para r que es mas convenientepara implementarse en una computadora es



17.1 REGRES IÓN L INEAL

E J E M P L O 17.2 Estimación de errores para el ajuste l ineal por mínimos cuadrados

Enunciado del problema. Calcule la desviación estándar total, el error estándar del

estimado y el coeficiente de correlación para los datos en el ejemplo 17.1.

Solución. Las sumatorias se realizan y se presentan en la tabla 17.1. La desviación

estándar es [véase ecuación (PT5.2)]

22 7143Í±1}ZÍ = 1.9457

7 - 1

y el error estándar del estimado es [véase ecuación (17.9)]

/2.9911 n „ „ ,s y / x = ^ T - J = 0.7735

Así, ya que sy/x < Sy, el modelo de regresión lineal tiene mérito. La mejora adicional se

puede cuantificar por [véase ecuación (17.10)]

22.7143 - 2.9911r" =

2 - — =0. 86 822.7143

r = V 0 . 8 6 8 = 0.932

Los resultados indican que el 86.8% de la incertidumbre original ha sido explicada por

el modelo lineal.

Antes de proceder con el programa de cómputo para regresión lineal, debemos to

mar en cuenta algunas consideraciones. Aunque los coeficientes de correlación propor

cionan una manera fácil para medir la bondad del ajuste, se deberá tener cuidado de no

darle más significado que el que ya tiene. Así como r es "cercana" a 1 no significa que el

ajuste es necesariamente "bueno". Por ejemplo, es posible obtener un valor relativamen

te alto de r cuando la relación en turno entre y y x no es lineal. Draper y Smith (1981)

proporcionan guías y material adicional con respecto al aseguramiento de los resultados

para regresión lineal. Además, como mínimo, usted debería siempre inspeccionar una

gráfica de los datos junto con su curva de regresión. Como se describe en la siguiente

sección, el software de métodos numéricos TOOLKIT incluye esas capacidades.

17 .1 .4 Progr ama de cómputo para regresió n lineal

Esto es una cuestión relativamente trivial para desarrollar un pseudocódigo para regre

sión lineal (véase figura 17.6). Como se mencionó antes, una opción de gráfica es crítica

para el uso efectivo e interpretación de regresión y se incluye en el software suplementa

rio de métodos numéricos TOOLKIT. Además, paquetes de software populares como

Excel y Mathcad pueden implementar regresión y tienen capacidades de graficación. Si

su lenguaje de computadora tiene capacidades de graficación, recomendamos que ex

panda su programa para incluir una gráf ica dey contra x mostrando ambos: los datos y la

línea de regresión. La inclusión de la capacidad resaltará mucho la utilidad del programu

en los contextos de solución de problemas.



4 7 4 REGRES I ÓN POR MÍ N I MOS CUAD RAD OS

5U 3 Regres(x, y, n, al, aO, syx, r2)

sumx = 0: sumxy = 0: st ~ 0

sumy = 0: sumx2 = 0: sr = 0

DO ¡ = 1n

sumx = sumx + x¡

sumy = sumy + y¡

sumxy = sumxy + x¡y¡

sumx2 = sumx2 + x^x,

EHDDO

xm = sumx/n

ym = sumy/n

a1 — (n*sumxy ~ sumx*sumy)/(n*sumx2 — sumx*sumx)

aO — ym — aUxm

DO i= 1,n

st = st + (y¡ - ymf

sr=sr+ (y¡ — a1*x¡ — aO) 2

EN D DO

syx = (sr/(n - 2)) 03

r2 = (st - s r j / s t

END Regres

F I G U R A 1 7 . 6

Algoritmo para regresión lineal.

EJ EM PL O 17 . 3 Reg res ión l inea l usa nd o la comp ut a dora

Enunciado de l problema. El paque te de so ftware de Métodos Num ér icos T OOL KITs adjunto a es te texto , contiene un program a de cómp uto para impl eme ntar regresión l i-: neal . Podem os usar es te sof tware para resolver un problem a de prueba h ipotético asocia-! do con la caída del paracaid is ta que se analizó en el capítu lo 1 . Un mo delo m atem áticoi teór ico para la velocidad del paracaid is ta fue dado como el s iguiente [véase ecuación

! (1 .10)] :

S " 1 (\ _ 0{-clm)t^C

{t) = — (1 - e

donde v = velocidad (m/s) , g = constante gravitacional (9 .8 m /s 2 ) , m = masa del paracaidista igual a 68.1 kg y c = coeficiente de arrastre de 12.5 kg/s. El modelo predice lavelocidad del paracaid is ta como una función del t iempo, com o se descr ibe en el ejem plo1.1. En el ejemplo 2.1 se desarrolló una gráfica de la variación de la velocidad.

Un modelo empír ico al ternativo para la velocidad del paracaid is ta es tá dado por

i;»; / l \



17.1 REGRESIÓN LINEAL

TABLA 1 7 . 2 Velocidades medidas y calculadas para la caída de l paracaidista,

v c a l c u l a d a v c a l c u l a d av m e d i d a . co n e l m o d e l o , co n e l m o d e l o ,

m / s m / s [ec . (1.10)] m / s [ ec . ( E l 7 . 3 . 1 ) ]T i e m p o , s o) b) c)

1 1 0 . 0 0 8 . 9 5 3 1 1 . 2 4 02 16 .30 1 6 . 4 0 5 1 8 . 5 7 03 2 3 . 0 0 2 2 . 6 0 7 2 3 . 7 2 94 2 7 . 5 0 2 7 . 7 6 9 2 7 . 5 . 5 65 3 1 . 0 0 3 2 . 0 6 5 3 0 . 5 0 9ó 3 5 . 6 0 3 5 . 6 4 1 3 2 . 8 5 57 3 9 . 0 0 3 8 . 6 1 7 3 4 . 7 6 68 4 1 . 5 0 4 1 . 0 9 5 3 6 . 3 5 1

9 4 2 . 9 0 4 3 . 1 5 6 3 7 . 6 8 710 4 5 . 0 0 4 4 . 8 7 2 3 8 . 8 2 911 4 6 . 0 0 4 6 . 3 0 1 3 9 . 8 1 612 4 5 . 5 0 4 7 . 4 9 0 4 0 . 6 7 813 4 6 . 0 0 4 8 . 4 7 9 4 1 . 4 3 714 4 9 . 0 0 4 9 . 3 0 3 4 2 . 1 1015 5 0 . 0 0 4 9 . 9 8 8 4 2 . 7 1 2

Suponga qu e a usted le gustar ía probar y com parar lo adec uado de esos dos m odelosmate mátic os. Esto se podr ía cumplir al me dir la velocidad real del paracaid is ta con valores conocidos de t iempo y comparar es tos resultados con las velocidades predichas deacuerdo con cada modelo .

Tal programa de colección de datos exper imentales se implemento , y los resultadosse enlis tan en la columna a) de la tabla 17 .2 . Las velocidades c alculadas para cada mod elo se enlis tan en las columnas b) y c).

Solución. La adecuidad de los mo delos se puede probar al graf icar la velocidad delmod elo calculad o contra la velocidad medid a. Se puede usar regresión l ineal para calcular la pendiente y el in tercepto de la gráf ica. Esta l ínea tendrá una pendiente de I , unin tercepto de 0 y una r 1 = 1 s i el modelo concuerda perfectamente con los datos . Unadesviació n significativa de esos valores se pue de usar com o un indicad or de la inadecuid adde l modelo .

La figura 17.7a y 17.76 son gráficas de la línea y datos para las regresiones de lasco lumnas b) y c) , respectivamente, contra la columna a) . Para el pr imer modelo [ecuación (1.10) com o se ilustra en la figura 17.7a]

« W wo = - 0 - 8 5 9 + 1.032u m e d l d a

y para el segundo modelo [ecuación (E17.3.1) como se ilustra en la figura 17.7/>]

"modelo = 5 - 7 7 6 + ° - 7 5 2 u m e d i d a

Esas grá f icas indican que la regresión l ineal entre los datos y cada uno de los modelos O N

altamente s ignif icativa. Ambos modelos ajustan los datos con un coef iciente de conciación mayor que 0 .99 .

I os modelos de prueba y selección son comunes y extremadiimenle imporliinlcNpara la realización de actividades en todos los campos de lu ingeniciiu, lil iniitcrinl pro




!

5 3 0

X

55

30 55

X

F I G U R A 1 7 . 7a) Resultados que usan regresión lineal para comparar las predicciones calculadas con elmodelo teórico [véase ecuación (1.10)] contra valores medidos, b) Resultados usandoregresión lineal para comparar predicciones calculadas con el modelo empírico [véaseecuación (El7.3.1] contra valores medidos.

porciona do en este capítu lo como antece dente, junt o con su sof tware, le proporcionará

una guía muy práctica para problemas de este t ipo .Sin embargo , el mod elo descr i to por la ecuación (1 .10) conforma para nuestra h ipó

tesis cr i ter ios de prueba mucho mejores que los descr i tos por la ecuación (E17.3 .1) yaque la pendiente y el in tercepto son m ás cercano s o casi igual a 1 y 0 . Así , aun que cadagráf ica está b ien descr i ta po ru ña l ínea recta, la ecuac ión (1 .10) parece ser mejor m odeloq u e l a ( E 1 7 . 3 . 1 ) .

Hay un defecto e on el anális is en el ejemplo 17.3 . El ejemplo no fue amb iguo, ynque el modelo emp ír ico [véase ecuación (E l 7 .3 .1) ] fue clara mente infer ior al de la ecuación (1 .10) . Así , la pendiente y el in tercepto para el model o emp ír ico fueron muc ho máscerc ano ! que el resultado deseado de 1 y 0, fue obvio cuál mo delo fue superior.



17 . 1 REGR ES I ÓN L I NEAL 477

Sin embargo, suponga que la pendiente fuera de 0.85 y que el intercepto lliorn de 2,Obvia mente esto har ía de la conclusión de que la pendiente y el in tercepto fueran I y 0,

un debate abier to . De m anera clara, más que recaer en un ju ic io subjetivo , ser ía prefer ib le basar tal conclusión sobre un cr i ter io cuanti tat ivo .

Esto se puede hacer al calcular los in tervalos de conf ianza para los parámetros delmodelo en la misma forma que desarro llamos los in tervalos de conf ianza para la mediaen la sección PT5.2 .3 . Regresaremos a es te punto al f inal del presente capítu lo .

1 7 . 1 . 5 L i n e a r i z a c i ó n d e r e l a c i o n e s n o l i n e a l e s

La regresión l ineal proporciona una técnica poderosa que ajusta a la "mejor" l ínea lo sdatos . Sin embargo, es tá predicha sobre el hecho de que la relación entre las var iablesdependientes e independientes es l ineal . Este no es s iempre el caso , y el pr imer paso encualquier anális is de regresión deber ía ser graf icar e inspeccionar en forma visual paraasegurarn os s i se puede usar un model o l ineal . Por ejemplo , la f igura 17 .8 mue stra a lgunos datos que son obviamente curvil íneos. En algunos casos, técnicas tales como regresión por polinomios, las cuales se descr iben en la sección 17.2 , son apropiadas. Paraot ros , se puede usar transformaciones para expresar los datos en una forma que se acompatib le con la regresión l ineal .

F I G U R A 1 7 . 8

o) Datos no adecuados para la regresión lineal por mínimos cuadrados, b) Indicación deque es preferible una parábola.

x

b)




Un ejemplo es el modelo exponencial

y = a i e h l X (1712)donde a x y b x son constantes . Este modelo se usa en muchos campos de la ingenier íapara caracter izar cantidades que aumentan (b { posit ivo) o d isminuyen (b ¡ negativo) uuna velocidad que es d irectamente proporcional a sus propias magnitudes. Por ejemplo ,el crecimiento poblacional o el decaimiento radiactivo pueden exhibir tal com por tamiento .Com o se i lustra en la f igura 17 .9a, la ecua ción representa una relación no l ineal (para /),

0) entre y y x.Otro ejemplo de modelo no l ineal es la s imple ecuación de potencias

y = a 2xhl (17.13)

F I G U R A 1 7 . 9o) La ecuación exponencial, b) la ecuación por potencias y c) la ecuación de razón de crecimiento saturado. Los incisos <

e] y f] son versiones linearizadas de estas ecuaciones producto de transformaciones simples.

y f

c)

1/y A

Pendiente • /y*j



17 .1 REGRE S I ÓN L I NEAL 47?

Enunciado del problema. Ajustar la ecuación (17.13) con los datos en la tabla 17.3mediante transformaciones logar í tmicas de los datos .

donde a 2 y b 2

son coeficientes constantes. Este modelo tiene amplia uplicnbiliclud 0 1 1todo s los cam pos de la ingeniería. Com o se ilustra en la figura 17.9/) , la ecu ació n (purab2^0o 1) es no lineal.

Un tercer ejemplo de un modelo no l ineal es la ecuación de razón de crecimientosaturado [ recuerde la ecuación (El7 .3 .1) ]

y = a

' l 7 ^ - <17

-14

>

donde a 3 y b 3 son coef icientes constantes . Este modelo , el cual es de manera par t icularmuy adecuado para caracter izar la razón de crecimiento poblacional bajo condicioneslimitadas, también representa una relación no l ineal entre y y x (véase f igura 17 .9c) queiguala o "satura", e n tanto x aumen ta .

Las técnicas de regresión no l ineal es tán d isponibles para ajustar esas ecuaciones a

datos exper imentales de manera d irecta. (Observe que analizaremos la regresión no l ineal en la sección 17.5 . ) Sin embargo, una al ternativa s imple es usar manipulacionesmatemáticas para transformar las ecuaciones en una forma l ineal . Después, se puedeemplear la regresión l ineal s imple para ajustar las ecuaciones a datos .

Por ejemplo , la ecuación (17 .12) se puede l inear izar al tomar su logar i tmo naturalpara dar

ln y = ln a\ + b\X ln e

Pero como ln e = 1,

ln y = ln ax

+ b\X (17 .15)

Así, una gráf ica de ln y contra x dará una línea recta con una pendiente de ¿>, y unin tercepto de ln a

l (véase la figura 17.9c/) .La ecua ción (17 .14) es linear izada al tomar su base logar i tmo 10 para dar

lo g y = b 2 lo g x + lo g a 2 (17 .16)

De este mo do, una gráf ica de log y contra log x dará una l ínea recta con una pendie nte deb

2y un in tercepto de log a

2 (figura 17.9e).La ecu ación (17 .14) es l inear izada al inver t ir la para dar

1 b 3 1 1- = + — (17.17)y « 3 X a 3

De esta forma, una gráfica de 1/y contra l/x será l ineal, con una pend iente de b-¡/a 3 y unin te r cep to de l /a 3 (véase la figura 17.9/) .

En sus contornos transformados, es tos modelos se ajustan mediante regresión l inealpara evaluar los coef icientes constantes . Podr ían ser de nuevo conver t idos en su estadoor ig inal y usados para propósitos predict ivos. El ejemplo 17.4 i lustra es te procedim ientopara la ecuación (17 .13) . Ade más , la sección 20.1 proporciona un ejem plo de ingen ier íade la misma clase de cálculo .

EJ EM PL O 17 . 4 L i nea r iza c i ón d e una ecua ción d e p otencia s



4 80 REGRESIÓN POR MÍNIMOS CUADRADOS

Solución. La figura 17. IOÍÜ es una gráfica de los datos originales en su estado no trans

formado . La figura 17.10¿> muestra la gráfica de los datos transfo rmados . Una reg resión

lineal de éstos mediante log dan el resultado

l o g v = 1.75 log x - 0.300

TABLA 17 .3 Datos que serán ajustados con

la ecuación de potencias.

x y log x log y

1 0.5 0 -0.301

2 1.7 0.301 0.2263 3.4 0.477 0.534

4 5.7 0.602 0.7535 8.4 0.699 0.922

l FIGURA 1 7 . 1 0

| a) Gráfica dedatos no transformados con la ecuación de potencias que ajusta los datos.

¡; b] Gráfica de datos transformados que se usan para determinar los coeficientes

de la ecuación de potencias.



17 . 2 REGRES I ÓN D E POL I NOMI OS 411

i Así , el in tercepto , log a 2, igual —0.300, y por tanto , al tomar el an ti logar i tmo, u¡ —

\ 1 0 _ 0 3 = 0 .5 . La pendiente es b 2 = 1 .75. En consecuencia, la ecuac ión de potencia s es

y • 0 . 5 . v l 7 S

Esta curva, como se gráf ica en la figura 17 .10a, indica un buen ajuste.

1 7 . 1 . 6 C o m e n t a r i o s g e n e r a l e s s o b r e r e g r e s i ó n l i n e a l

Antes de proceder con regresión curvil ínea y l ineal múlt ip le, deb em os enfatizar la naturaleza in troductor ia del mater ial an ter ior sobre regresión l ineal . No s hem os conc entrado

en la der ivación s imple y uso práctico de ecuaciones para ajustar datos . Deber ía es tarconsciente del hecho de que hay aspectos teór icos de regresión que son de impor tanciapráctica, pero que van más al lá del alcance de este l ibro . Por ejemplo , alguna s su posiciones estadís t icas que son inherentes en los procedimientos por mínimos cuadrados l ineales son

1 . C ada x tiene un valor fijo; no es aleatorio y es conocido sin error.2. Los valores y son var iables aleatorias independientes y todas t ienen la mism a var ianza.3. Los valores de y para una x dada deben ser normalmente d is tr ibuidos.

Tales suposiciones son relevantes para la der ivación adecuada y uso de regresión .Por ejemplo , la pr imera suposición s ignif ica que 1) los valores x deben estar l ibres deerrores y 2) la regresión de y contra x no es la misma que la de x contra y (pruebe el

problema 17.4 al f inal del capítu lo) . Usted debe consultar o tras referencias tales comoDraper y Smith (1981) para apreciar aspectos y matices de regresión que están más al ládel alcance de este l ibro .

1 7 . 2 R E G R E S I Ó N D E P O L I N O M I O S

En la sección 17.1 se desarro lló un procedimie nto para obtener la ecuación de una l ínearecta por me dio del cr i ter io de mínim os cuadrad os. Algu nos datos de ingenier ía , aunqueexhiben un patrón ma rcad o com o el que se v io en la f igura 17 .8 , es tá pobrem ente re presentado por una l ínea recta. Para esos casos, una curva podr ía ser más adecu ada para elajuste de los datos . Com o se analizó en la sección anter ior , un mé todo para cum plir co neste objet ivo es usar transformacion es. O tras al ternativas son ajustar polinom ios con losda to s med ian te regresión de polinomios.

El procedimiento de mínimos cuadrados se puede fácilmente extender al ajuste dedatos con un polinomio de orden super ior . Por ejemplo , suponga que ajustamos unpolinomio de segundo orden o cuadrático:

y = a 0 + a\X + atx 1 + e

Para este caso la suma de los cuadrados de los res iduos es [compare con la ecuación(17.3)]

Sr = ]T] (y i ~ «o - a\x¡ - a

2xf)

(17 .18)



4 82 REGRESIÓN PO R MÍNIMOS CUADRADOS

Sigu iendo el p r o c e d i m i e n t o de la secc ión an te r io r , tomamos la der ivada de la ecuación (17 .18) con respecto de cada uno de los coef icientes desconocidos del polinomio,c o m o en

Estas ecuaciones se pueden igualar a cero y reordenar para desarro llar el siguiente con

j u n t o de ecuac iones no rmal :

donde todas las sumator ia s son de i = 1 hasta n. Observe que las tres ecuaciones a nter iores son l ineales y t ienen tres incógnitas : aQ, ai y a2. Los coef icientes de las incógnitas se

pueden evaluar de manera d irecta a par t ir de los datos observados.Para este caso , vemos que el prob lema de dete rminar un po l inomio por m í n i m o s

cuad rados de segundo orden es equivalente a resolver un s i s tema de t res ecuacioneslineales s imultán eas . Las técnicas para resolver tales ecuacione s fueron analizada s en la

par te tres .El caso en dos d imensiones puede extenderse con facil idad a un po l inomio de m-

és imo o rden como

y = a0 + a\x + a2x2 H h amx

m + e

El anális is anter ior se puede fácilmente extender a este caso más general . Así , podemosreconoc er que la determ inación de los coef icientes de un polino mio d e /n-ésimo orden esequivalente a resolver un sis tema de m + 1 ecuaciones l ineales s imultáneas . Para es tecaso , el error es tándar se formula como

Esta cantidad es dividida entre n — (m + 1) , ya que (m + 1) coef icientes obtenidos de

los datos (a0, a , , . . . , am) se usaron para calcular S,.; así , hemos perdido m + 1 grados de

l iber tad . Ad em ás del er ror es tándar , un coef iciente de determinación puede ser calculado para una regresión de polinomios con la ecuación (17 .10) .

(17 .19)

(17.20)

E J E M P L O 1 7 . 5 Regresión de polinomios

E n u n c i a d o del prob lema. A jus ta r a un po l inomio de segundo orden los datos en I U N

dos primeras columnas de la tabla 17 .4 .



17.2 R E G R E 5 I O I N DE P O L I N O M I O S

TABLA 17.4 Cálculos para un análisis de error del ajustecuadrático por mínimos cuadrados.

x ¡ y , ( y ; - y ) 2 ( y / - O O - O I X f - o 2 x ? )

0 2 .1 5 4 4 . 4 4 0 . 1 4 3 3 2

1 7.7 3 1 4 . 4 7 1 . 0 0 2 8 62 1 3 . 6 1 4 0 . 0 3 1 . 0 8 1 5 83 2 7 . 2 3 . 1 2 0 . 8 0 4 9 14 4 0 . 9 2 3 9 . 2 2 0 . 6 1 9 5 15 6 1 . 1 1 2 7 2 . 1 1 0 . 0 9 4 3 9

2 1 5 2 . 6 2 5 1 3 . 3 9 3 . 7 4 6 5 7

F IGURA 17.11

Ajuste de un polinomio de segundo orden.

Solución. A partir de los datos dados,

m = 2

X> =

5

X>=

79

n = 6

í>-152.6 = 585.6

X = 2.5 2>?-55 = 2 4 8 8 .

y = 25.433 225

Por tan to , las ecuaciones lineales simultáneas son

" 6 15 55 " « 0 152.6

15 55 225 0\ • = • 585.6

55 225 979 «2 2488.8



4 8 4 REGRES I ÓN POR MÍ N I MOS CUAD RAD OS

Resolviendo estas ecuaciones con una técnica tal como la el iminación de Gauss se t ienea 0 = 2.47857, a , = 2 .35929 y a 2 = 1 .86071. Por tanto , la ecua ción cuadrática pormínimos cuadrados para es te caso es

y = 2 .47857 + 2 .35929* + 1 .86071* 2

El er ror es tándar del es t imado con base en la regresión de polinom ios es [véase ecuación(17.20)]

Sy/,3.74657

= 1.126 - 3

El coef iciente de determinación es

, 2 5 1 3 . 3 9 - 3 .7 46 5 7r2 = = 0 .99851

2513 .39y el coef iciente de correlación es r — 0 .99925 .

Estos resultados indican que el 9 9 . 8 5 1 % de la incer t idumbre or ig inal la resolv ió elmod elo . Este resultado sopor ta la conclusión de que la ecuación cuadrátic a representa unexcelente ajuste, como es también evidente de la f igura 17.11.

1 7 . 2 . 1 A l g o r i t m o p a r a r e g r e s i ó n d e p o l i n o m i o s

Un a lgor i tmo para regresión de polinom ios es expuesto en la f igura 17 .12. Observe quela pr incipal tarea es la generación de los coef icientes de las ecuaciones normal [véaseecuación (17 .19)] . (El pseudocódig o para el cumplim iento de esto se halla presente en la

f igura 17 .13.) Entonc es, las técnicas de la par te tres se pueden a plicar para resolver es tasecuaciones s imultáneas para los coef icientes .Un problema potencial asociado con la implementación de regresión de polinomios

en la computadora es que las ecuaciones normales algunas veces están mal condicionadas . Esto es en par t icular cier to para vers iones de orden super ior . Para esos casos, loscoef icientes calculados podr ían ser al tamente susceptib les al er ror de redondeo y , enconsecuencia, los resultados pueden ser inexactos . Entre o tras cosas , es te problema serelaciona con la es tructura de las ecuaciones normal y por el hecho de que para lospolinomios de orden super ior las ecuaciones normales pueden tener coef icientes muygrandes y mu y peque ños. Esto se debe a los coef icientes y sum ator ias de los datos elevados a potencias .

F I G U R A 1 7 . 1 2Algoritmo para la implementación de polinomios y regresión lineal múltiple.

P a s o 1 : Introduzca el orden del polinomio sujeto a ajuste, m.

P a s o 2 : Integre el número de datos, n.

P a s o 3 : Si n < m + 1, imprima un mensaje de error que indique que la regresión no es posible

y termine el proceso. Si n S m + 1, continúe.P a s o 4 : Calcule los elementos de la ecuación normal en la forma de una matriz aumentada.

P a s o 5 : Resuelva la matriz aumentada para los coeficientes a 0 , d i , a 2 , . . ., a m , por m edio de un

método de eliminación.P a s o d i Imprima los coeficientes.



17 . 2 REGR ES I ÓN D E POL I NOMI OS 4 1 9

DO i = 1, o rd í r + 1PO = 1, i

eum — 0

DOi = 1, n

eum — eum + ¿£

EN D DO

a¡j = eum

ap — eum

EN D DO

eum = 0

D0£ = 1,n

eum — eum + y ( • x£ ~ 1

EN D DO

ai, arder + 2 = 5 U M

EN D DO

F I G U R A 1 7 . 1 3Pseudocódigo para ensamblar los elementos de las ecuaciones normal para regresión depolinomios.

Aunque las estrategias para disminuir el error de redondeo analizado en la partetres, tal como el pivoteo, pueden ayudar a remediar en forma parcial este problema, unaalternativa más simple es usar una computadora cotí más aira precisión. Por fortuna, la

mayoría de los problem as prácticos están limitados a polino mi os de orden inferior paralos cua les el error de redondeo es insignificante. En situaciones d onde se requieren versiones de orden superior, se dispone de otras alternativas para ciertos tipos de datos. Sinembargo, esas técnicas (tal como polinomios ortogonales) están más allá del alcance deeste libro. El lector debería consultar textos sobre regresión, com o el d e Draper y Smith(1981), para información adicional con respecto al problema y posibles alternativas.

EJ EM PL O 1 7 . 6 Reg res ión de p o l i nomi os p or med i o d e la comp ut a dora

Enunciado del problema. En el sof tware de Métodos Nu mérico s TOOLKIT adjunto aeste libro se t iene un programa de computadora de uso amigable para implementar laregresión de polinomios. Se puede usar este software para el ajuste de polinomios conlos siguientes datos:

X 2 4 5 6 6 7 9 1 0.5 7. 5

y 6 2 3 7 8 8 1 5 3 7

Solu ción . Presione el ajuste de datos con el botón Curve sobre el men ú principal delTOO LKIT para obtener una pantalla en blanco similar a la de la f igura 17.14. Esta pan-

I talla contiene espa cios para la entrada y salida de información n ecesaria para ajustar lo sI datos con un polinomio de regresión por mínimos cuadrados de nj-ésimo orden,J El primer pa so es presionar lo s valor es de entrada X contra Y en la tabla e introducirI hasta 100 pares de valores para X y Y. De spué s usted podría decidir graficar los da tos



4 84 REGRES I ÓN POR MÍ N I MOS CUAD RAD OS

Resolviendo estas ecuaciones con una técnica tal como la el iminación de Gauss se t ienea 0 = 2 .47857, a x = 2 .35929 y a 2 — 1.86071. Por tanto , la ecuación cuadrática pormínimos cuadrados para es te caso es

y = 2 .47857 + 2 .35929x + 1 .86071* 2

El er ror es tándar del es t imado con base en la regresión de polinom ios es [véase ecuació n(17.20)]

3 .74657= 1.12

6 - 3El coef iciente de determinación es

, 2513.39 - 3 .74657r2 = = 0 .99851

2513 .39y el coef iciente de correlación es r = 0 .99925 .

Estos resultados indican que el 9 9 . 8 5 1 % de la incer t idumbre or ig inal la resolv ió elmod elo . Este resultad o sopor ta la conclusión de que la ecuació n cuadrática representa unexcelente ajuste, como es también evidente de la f igura 17.11.

1 7 . 2 . 1 A l g o r i t m o p a r a r e g r e s i ó n d e p o l i n o m i o s

Un algor i tm o para regresión de polino mios es expuesto en la f igura 17 .12. Observ e quela pr incipal tarea es la generación de los coef icientes de las ecuaciones normal [véaseecuac ión (17 . 19) ] . (El pseudocódigo pa ra el cumplim iento de esto se halla presente en la

f igura 17 .13.) Entonc es, las técnicas de la par te tres se pued en aplicar para resolver es tasecuaciones s imultáneas para los coef icientes .

Un problema potencial asociado con la implementación de regresión de polinomiosen la computadora es que las ecuaciones normales algunas veces están mal condicionadas . Esto es en par t icular cier to para vers iones de orden super ior . Para esos casos, loscoef icientes calculados podr ían ser al tamente susceptib les al er ror de redondeo y , enconsecuencia, los resultados pueden ser inexactos . Entre o tras cosas , es te problema serelaciona con la es tructura de las ecuaciones normal y por el hecho de que para lospolinomios de orden super ior las ecuaciones normales pueden tener coef icientes muygrandes y mu y peque ños. Esto se debe a los coef icientes y sumator ias de los datos elevados a potencias .

F I G U R A 1 7 . 1 2Algoritmo para la implementación de polinomios y regresión lineal múltiple.

P a s o 1 : Introduzca el orden del polinomio sujeto a ajuste, m .P a s o 2 : Integre el número de d atos, n.P a s o 3 : Si n < m + 1, imprima un mensaje de error que indique que la regresión no es po sible

y termine el proceso. Si n 5 : m + 1, continúe.P a s o 4 : Calcule los elementos de la ecuación normal en la forma de una matriz aumenlada.

P a s o 5 : Resuelva la matriz aumentada para los coeficientes o 0 , a i, a 7 , . . ., a,„, por medio do un

método de eliminación.

P A S O 6 L Imprjma los coeficientes.



17 . 2 REGRES I ÓN D E POL I NOMI OS 4 8 3

DO i = 1, order + 1

DOj = 1, i

k=i + j - 2

eum = 0

DOt =\n

eum = eum + >¿€

EN D DO

a, • = eum

ap = eum

EN D DO

eum = 0

D 0€ = 1,n

eum — eum + y(

• x ¿-

'EN D DO

ai. arder + 2= 5U™EN D DO

F I G U R A 1 7 . 1 3Pseudocódigo para ensamblar los elementos de las ecuaciones normal para regresión depolinomios.

Aunque las estrategias para disminuir el error de redondeo analizado en la partetres, tal como el pivoteo , pueden ayudar a remediar en forma parcial este problema, unaalternativa m ás simple es usar una computadora co r/ más alta precisión . Por fortuna, la

mayoría de los problem as prácticos están limitad os a polinprriios de orden inferior paralo s cuales el error de redondeo es insignificante. En situaciones donde se requieren versione s de orden superior, se dispone de otras alternativas para ciertos t ipos de datos. Sinembargo, esas técnicas (tal como polinomios ortogonales) están más allá del alcance deeste libro. El lector debería consultar textos sobre regresión, com o el de Draper y Smith(1981), para información adicional con respecto al problema y posibles alternativas.

EJ EM PL O 17 . 6 Reg res ión d e p o l i nomi os p or med i o d e la comp ut ad ora

Enunciado del problema. En el software de Métod os Num érico s TOOLKIT adjunto aeste libro se t iene un programa de computadora de uso amigable para implementar laregresión de polinomios. Se puede usar este software para el ajuste de polinomios conlos siguientes datos:

2 4 5 6 6 7 9 1 0.5 7.5

Solu ción . Presione el ajuste de datos con el botón Curve sobre el men ú principal delTOOL KIT para obtener una pantalla en blanco similar a la de la f igura 17.14. Esta pantalla contiene espacios para la entrada y salida de información necesaria para ajustar losdatos con un polinomio de regresión por mínimos cuadrados de m-ésimo orden.

El primer paso es presionar los valores de entrada X contra Y en la tabla e introducirhasta 10 0 pares de valores para X y Y . Después usted podría decidir granear los datos



486 RKMMIÓN PO R MÍNIMOS CUADRADOS

solos antas de realizar decis iones con respecto al orden del polinomio. Esto se hacemed iante un proced imien to s imilar al descr i to en el ejemplo 2 .1 . La inspecc ión de losdatos muestra dos p icos y sugiere que un polinomio de al menos cuar to orden ser ía eladecuado. Para nuestro ejemplo , pr imero in tentaremos un polinomio de quinto orden.Simplemente in troduzca un valor de 5 para el orden del polinomio y graf ique los parámetros en el cuadro Entrada de parámetros y haga cl ic en los botones de red Cale y Plo t(en el proceso c amb ian los botones a un colo r negro) para producir la f igura 17 .14. Laforma para determinar el mejor orden se puede explorar al examinar cóm o el er ror es tándarvar ía como una función del orden de regresión . Los resultados para var ios órdenes deregresión en el ajuste se tabula en la s iguiente pág ina:

FIGURA 1 7 . 1 4Pantalla del TOOLKIT de métodos numéricos para una regresión polinomial de quinto orden.

P^iwaete* ~ Valué:| O Í d e t of Paíy h

í Plot Xmtn 0Delta X 1Plot Ymin 0Delta Y 1

l i i l l l i +

'i ' 12 6

3 2 4 2, 3 5 3

4 6 7

5 6 86 7 8, 7 9 1

8 1 59 .5 3

V 10 7. 5

I ;.>,;.«8SftÉ!£ ; ¿J.:.

X - 3 * > fiettift ,.. . . Wm...:,Y « 2.304472 Standard Erior 1 000334

Coef of Deter .3332087Coir Coef .3660273

Oth oider coef -6.573234

J B É B J ' C S 3 Q 5 E 3 * ' C H C 3 - i m« i XMjml

FIGURA 1 7 . 1 5Gráfica de una regresión polinomial de octavo orden.

Vafc»Ordet of Poly 8

Plot Xmro 0Delta X 1 1Plot Ymin 0Delta Y 1 i

l w p u t X v * ¥ V a f o e *

x v1 2 6

i 2 4 2

í 3 5 3

1 4 6 7í 5 6 8

1 6 7 . 8Ü 7 9 1

8 1 59 .5 3

10 75 7

CafcYÍOT InputX

ftewft Vafe» 'Standatd Enot 1 154277Coef of Detei .9777941

Con Coef .9888347Oth order coef -1 057965 | £



17.3 REG RESIÓ N L INEAL MÚLT IPLE 4 87

Orden 2 3

A6 /

Error estándar 2.7 1 2 .6 9 2 .3 4 1.38 1 . 0 0 1.12 1 . 1 7

Obse rve que el er ror es tándar cae de manera dramática d el orden 3 al 4 y alcanza unmínimo para el polinomio de orden 5 . Esto sugiere que no se ha ganado mucho al gastaren esfuerzo de cómputo para ejecutar la regresión más al ta que la de quin to orden.

La f igura 17 .15 muestra las gráf icas para el caso de un octavo orden. Para este c aso ,lo s pun tos ex t r emo empiezan a se r un p rob lema en una manera s imi la r a la de lain terpolación de orden super ior (analizarem os este fenómeno c on más detal le en el próx imo capítu lo) . La f igura 17 .15 mue stra que el polinom io de octavo orden produce valoresde Y negativos para valores de X entre 8 y 9 . Observe tam bién en las f iguras 17 .14 y17.15 que mientras las curvas de regresión s iguen la tenden cia de los datos , es al tam enteinapropiado extrapolar los valores Y má s al lá del rango de los datos para X .

La in terpo lación se puede ejecutar al in troducir un valor en la tabla para X en el C aleY para una X in troducida. Por ejemplo , en X = 3, Y = 2 .304472 como calculada con elpolinom io de quinto orden (véase la f igura 17 .14) .

Por ú lt imo, demos un v is tazo a la tabla de resultados en la par te derecha infer ior .Los pr imeros tres resultados son resúmenes estadís t icos de la regresión: er ror es tándar ,coef iciente de determinación y coef iciente de correlación . Observe cómo esos valorescambian para d iferentes órdenes de regresión . La barra de despliegue sobre la tabla deresultados se usa para observar los coef icientes reales de la regresión de polinom ios. D enuevo, esos valores cambian con d iferentes órdenes.

Una extensión ú ti l de la regresión l ineal es el caso donde y es una función lineal de dos omás var iables independientes . Por ejemplo , y podría ser una función lineal de x

l y x2,

como en

y = c ? o + « i A ' i + 2X2 + eTal ecuació n es en par t icular ú t il cuando se ajustan datos exper im entales donde la var iab le su jeta a es tud io es a me nudo una función de o tras dos var iables . Para es te caso en do sdim ension es, ta " l ínea " de regresión pasa a ser un "plan o" (véase la figura 17 .16) .

Com o en los casos anter iores , los "mejo res" valores de los coef icientes son dete rmi

nados al realizar la suma de los cuadrados de los res iduos,

R E G R E S I Ó N L I N E A L M Ú L T I P L E

(17 .21)

y d iferenciando c on respecto a cada uno de los coef icientes descon ocidos,



488 R E G R E S IÓ N P O R M Í N I M O S C U A D R A D O S

F I G U R A 1 7 . 1 6

Ilustración gráfica de regresión lineal múltiple donde yes una función lineal de x, y x 2 .

dS r

3 (3 ?= - 2 ^ x 2 / (y, - í/ 0 - a i x a ~ « 2 * 2 / )

Los coef icientes dan la suma mínima de los cuadrados de los res iduos y se obtienen aligual las der ivadas parciales a cero y expresando el resultado en forma de matr iz como

n

E x , /E x 2 í Exi, -x 2l

E x 2 /

E x 2,. T,xi¡x2¡

E x 2

a 0 Ey,-ai = • Exi,y,-a 2 S x 2 í y , '

(17 .22)

¡ EJEMPLO 1 7 . 7 Reg res ión l inea l mú l ti p le

Enunciado del problema. Los s iguientes datos se calcularon con la ecuación y = 5 +

4x , — 3 x 2 :

x l X2y

0 0 52 1 102.5 2 91 3 04 6 37 2 2 7

Use regresión lineal múltiple para ajustar estos datos.



17.3 REGRESIÓN LINEAL MÚLTIPLE

Solución. Las sumatorias requeridas para desarroll ar la ecuación (17.22) se calculanen la tabla 17.5. El resultado es

6 16.5 14

16.5 76.25 48

14 48 54

la cual se puede resolver mediante un método como el de eliminación de Gauss para

ao = 5 ai = 4 a2 = — 3

que es consistente con la ecuación original a partir de la cual los datos se derivaron.

TABLA 17.5 Cálculos requeridos para desarrollar las ecuaciones normal para

el ejemplo 17.7.

y * i X

2X 2

A x , x 2 * i y x 2 y

5 0 0 0 0 0 0 0

10 2 1 4 i 2 20 109 2.5 2 6.25 4 5 22.5 18

0 1 3 1 9 3 0 0

3 4 6 16 36 24 12 18

27 7 2 49 4 14 189 54

1 54

16.5 14 76.25 54 48 243.5 100

a0 54

a\ • = • 243.5

a2100

El caso anterior en dos dimensiones se puede fácilmente extender a m dimensiones,

como en

y = an + a\X\ + a2x2 -| Y amxm + e

donde el error estándar se formula como

S y , X ~\¡n-(m + l)

y el coeficiente de determinación se calcula como en la ecuación (17.10). En la figura

17.17 se enlista un algoritmo para preparar las ecuaciones normal.

Aunque hay ciertos casos donde una variable está linealmente relacionada con otras

dos o más variables, la regresión múltiple tiene utilidad adicional en la derivación de

ecuaciones de potencias de la forma general

y=a 0x 1

1x 2

2---xm"

Tales ecuaciones son extremadamente útiles cuando se ajustan datos experimentales.

Para usar regresión lineal múltiple, la ecuación se transforma al tomar su logaritmo para

obtener

log y = log üo + a\ log x\ + a2 log xj-i Y am log xm



REGRES I ÓN POR MÍ N I MOS CUAD RAD OS

DO I = 1, order + 1

DOj = 1, i

eum — 0

DOÍ= 1,n

eum = eum + x¡-^( • XJ-K

END DO

a u = eu m

a- s — eu m

END DO

eum — 0

DOl = \,n

eum = eum + y e • x¡_ 1(

END DO

A¡,CRDER+2 = 5 U M

END DO

F I G U R A 1 7 . 1 7Pseudocódigo para ensamblar los elementos de las ecuaciones normal para regresiónmúltiple. Observe que además de guardar las variables independientes en X ] , x 2 ¡, etcétera,los 1 se deben guardar en X Q ,• para trabajar este algoritmo.

Esta transform ación, es similar en ese ncia a la que se usó e n la secció n 17.1.5 y en elejemplo 17.4 para ajustar a una ecuación po r potencia s cua nd oy fue una función de unasola var iable x. La sección 20.4 proporciona un ejem plo de tal ap licación p ara dos var iab les independientes .

1 7 . 4 F O R M A G E N E R A L L I N E A L P O R M Í N I M O S C U A D R A D O S

Hasta es te punto nos hemos concentrado en la mecánica de obtención de ajustes pormín imo s cuadrados para algunas funciones s im ples con datos . Ante s de cambiar a regres ión no l ineal , hay var ios puntos que nos gustar ía analizar para enr iqu ecer nuestra comprensión del mater ial precedente.

1 7 . 4 . 1 F o r m u l a c i ó n g e n e r a l d e u n a m a t r i z p a r a m í n i m o s

c u a d r a d o s l i n e a l e s

En págin as anter iores hem os in troducido tres t ipos de regresión: l ineal s imple, polinomialy l ineal múlt ip le. De hecho, es tas tres per tenecen al s iguiente modelo general de mínimos cuadrados l ineales:

y — a r j z o +a \ z \ + a 2 z 2 ^ 1" amzm + e (17 .23)

donde z 0 , z , , . . . , zm son las m + 1 funciones d iferentes . Se pued e ver con facilidad cómola regresión l ineal s imple y múlt ip le encajan dentro de este modelo ; es decir z () — 1, z,j C | , z 2 = x 2, • •., z„, = xm. Además, la regresión de polinomios se incluye también si Iris ;son monomios simples como en z 0 = xü — I , z, = x, Z 2 = x 2,. .., Z M — x" 1.



Observe que la terminología " l ineal" se ref iere só lo a la dependencia del modelosobre sus parámetros (es decir , las a) . Como en el caso de regresión de polinomios, lasmismas funciones pueden ser al tamente no l ineales . Por ejemplo , las z pueden sers inusoidales , como en

y = a0 + a, eos (coi) +a2 se n (coi)Tal formato es la base del anális is de Four ier descr i to en el capítu lo 19 .

Por o tro lado, un modelo de apar iencia s imple como

f(x) = a 0 (1 - e-"'*)

es cier tamente no l ineal porque no puede ser manejado en el formato de la ecuación(17.23) . Regresaremos a tales modelos al f inal de este capítu lo .

Mientras tanto , la ecuación (17 .23) se puede expresar en notación matr icial como

{Y} = [Z]{A] + {E} (17 .24)

donde [2] es una m atr iz de los valores calculados de las funciones z en los valores me didos de las var iables independientes ,

[Z] =

Z 0 2 Z l2

Zrj« z l «

zm 1

zm2

donde m es número de var iables en el mo delo y n es el número de datos . Como n>m +

1, usted deber ía rec onoce r que la mayor ía de las veces [Z] no es una matr iz cuadrada .El vector columna {Y} contiene los valores observados de la var iable dependiente

{Y}T = [yi yi ••• y„JEl vector columna {A } contiene los coef icientes desconocidos

{A } T — [a 0 a¡ ••• a m \

y el vector columna {E } contiene los res iduos

[E}T = \_e\ e2 ••• e„ \

Com o se realizó a través de este capítu lo , la suma de los cuadra dos de los re s iduospara este modelo se pueden def in ir como

n / m \ 2

Sr = Y,\x-Y,a¡zAi = \ \ j=0 I

Esta can tidad se puede m inimizar al tomar su der ivada parcial con respecto a cada uno d elos coeficientes y fijar los resultados de la ecuación igual a cero. La salida de este proceso son las ecuaciones normal que se pueden expresar brevemente en forma de matr izc o m o



4 9 2 REGRESIÓN PO R MÍNIMOS CUADRADOS

[lZ]r

[Z]]{A} = {[Z]T

{Y}}Se puede demostrar que la ecuación (17 .25) es , de hecho, equivalente a las ecuacionesnormal desarro lladas antes para regresión l ineal s imple, polinomial y múlt ip le.

Nuestra pr incipal motivación para las anter iores ha s ido i lustrar la unión entre lostres procedimientos y mostrar cómo se pueden expresar de manera s imple en la mismanotación matr icial . También acondiciona la etapa para la s iguiente sección donde asimilaremos algo de conocimiento en las es trategias prefer idas para resolver la ecuación(17.25) . La notación matr icial tendrá también relevancia cuando veamos regresión nolineal en la ú l t ima secc ión de este ca pítu lo .

1 7 . 4 . 2 T é c n ic a s d e s o l u c i ó n

En los anális is anter iores en este capítu lo hemos encubier to el tema de las técnicas numér icas específ icas para resolver las ecuaciones normales . Ahora que hemos establecidola unión entre los d iversos modelos , podemos explorar es ta cuestión con mayor detal le .

Pr im ero , deber ía quedar claro que Gauss-Se idel no pue de usarse aquí debido a que lasecuaciones normal no son d iagonalmente dominantes . De esta manera dejamos a un ladolos métodos de el im inación. Para los actuales propósitos , podem os d iv id ir esas técnicas entres categor ías: 1) métodos de descomposición LU , inc luyendo e l iminac ión de Gauss ,2) método de Cholesky y 3) procedimiento de inversión de matr ices . Obviamente haytraslapes en esta clas if icación. Por ejemplo , el méto do de Choles ky es , de hecho, una descompos ic ión L U , y todos los procedim ientos se pueden formula r de tal forma qu e puedengenerar la matr iz inversa. Sin embargo, es ta clas if icación t ien e su mér ito en cada categor íay ofrece bene f icios con respecto a la so lución de las ecuaciones n orm ales .

Descompos ic ión LU. Si usted está in teresado sólo en aplicar un ajuste por mínimoscuadrados para el caso donde se conoce de antemano el modelo adecuado, cualquiera delo s p roced imien tos de descompos ic ión LU descr i tos en el capítu lo 9 son perfectamenteaceptables . De hecho, se puede también emplear la formulación de una descomposiciónno Lt/de el iminación de Gauss . Ésta es una tarea de programación relat ivamente d irectapara incorporar cualquiera de estos procedimientos en un algor i tmo de mínimos cuadrados l ineales . De he cho, s i se ha seguido un en foque m odular , es to es casi t r iv ial .

M éto do de Cholesky. El algor i tm o de descom posición de Choles ky t iene var ias ventajas con respecto a la so lución del problema general de regresión l ineal . Pr imero , es táexpresamente d iseñado para resolver matr ices s imétr icas como las ecuaciones normal.

De este mod o, es rápido y requiere menos espacio de almac enam iento para resolver taless is temas. Segundo, es idealmente adecuado para casos donde el orden del modelo [esdecir, el valor de m en la ecuación (17 .23)] no es conocido de an tema no (véase Rals ton yRabinowitz, 1978) . Un caso sujeto a tratamiento ser ía la regresión de polinomios. Paraeste caso , podr íamos saber a priori si un polinom io l ineal cuadrático , cúbico o de ordensuper ior es el "me jor" m odelo para descr ib ir nuestros datos . Debido tanto a la forma enla cual las ecuaciones normales se construyen como a la manera en la que procede elalgor i tmo de Cholesky (véase f igura 11 .3) , podem os desarro llar en forma sucesiva modelos de orden super ior de un modo en extremo ef iciente. En cada paso podr íamos examinar la suma residual del error de los cuadrados (¡y una gráfica!) para examinar si lainclusión de términos de orden super ior mejoran de ma nera s ignif icativa el ajuste.



i fA F O R M A GENERAL LINEAL POR MÍNIMOS CUADRADOS

493

La situación análoga para regresión lineal múltiple ocurre cuando se agregan variu-bles independientes, una a la vez, al modelo. Suponga que la variable dependiente de

interés es una función de un número de variables independientes: por ejemplo, tempera

tura, contenido de humedad, presión, etcétera. Podríamos primero realizar una regresión

lineal con la temperatura y calcular un error residual. En seguida se podría incluir el

contenido de humedad para realizar una regresión múltiple de dos variables y ver si la

variable adicional resulta en una mejora al ajuste. El método de Choles ky hace eficiente

el proceso, ya que la descomposición del modelo lineal podría solamente ser añadido al

incorporar una nueva variable.

Procedimiento de la matriz inversa. De la ecuación (PT3.6), recuerde que la matriz

inversa se puede emplear para resolver la ecuación (17.25), como en

{A} = [[Z]T

[Z]YL {[Z]T{Y}} (17.

Cada uno de los métodos de eliminación se puede usar para determinar la inversa y, así,

pueden ser usados para implementar la ecuación (17.26). Sin embargo, como aprendi

mos en la parte tres, éste es un enfoque ineficiente para resolver un conjunto de ecuacio

nes simultáneas. Así, si estuviéramos justamen te interesados en resolver los coefientes

de regresión, es preferible utilizar la aproximación de descomposición LU sin inversión.

No obstante desde una perspectiva estadística, hay un número de razones por las cuales

podríamos estar interesados en obtener la inversa y examinar sus coeficientes. Esas ra

zones se analizarán después.

17.4.3 Aspectos estadísticos de la teoría de mínimos cuadradosEn la sección P T 5 . 2 . 1 , revisamos un número de estadística descriptiva que puede usarse

para describir una muestra. Aquéllas incluyen la media aritmética, la desviación estándar

y la varianza.

Además de obtener una solución para los coeficientes de regresión, la formulación

de la matriz de la ecuación (17.26) proporciona estimaciones de su estadística. Se puede

demostrar (Draper y Smith, 1981) que la diagonal y los términos fuera de la diagonal de

la matriz [[Z]T [ Z ] ] _ 1 dan, respectivamente, las varianzas y las covarianzas 1 de las a. Si

los elementos de la diagonal de [[Z]r [Z]]~1 son designados como z j l , entonces

var (ai-i) = zj¡l52

y/x (17-27)

ycov (ai-i,aj)=zr}ijs*/x ( 1 7 2 g )

Estas estadísticas poseen un número de aplicaciones importantes. Para nuestros ac

tuales propósitos, ilustraremos cómo se pueden usar para desarrollar intervalos de con

fianza para el intercepto y la pendiente.

'Lacovarianza es una estadística que mide la dependencia de una variable con ot ra . Asi, cov (x,y) indica

la dependenciaácxyy. Por ejemplo, cov(x,y) 0 podría indicar queje yy son totalmente indcpcn



494 REGRESIÓN POR MÍNIMOS CUADRADOS

Mediante un enfoque similar al visto en la sección l'T'5.2.3, se puede demostrar que

los límites inferior y superior del intercepto se pueden formular como (véase Milton y

Arnold 1995 para más detalles)

<•><) — ta/2.,¡-2í(«o) U = fln + í a / 2 . H - 2 - V ( « ( ) ) (17.29)

donde s(a¡) — error estándar del coeficiente a- = K vai(ap. De manera similar, los lími

tes inferior y superior de la pendiente se pueden formular como

L — a\ - ta/2.„-is(a\) U = a\+tai2.n-2s(ci\) (17.30)

El siguiente ejemplo ilustra cómo esos intervalos se pueden usar para hacer inferencias

cuantitativas relacionadas con la regresión lineal.

EJ EM PL O 17 .8 Intervalos de confianza para regresión lineal

\ Enunciado del problema. En el ejemplo 17.3 usamos regresión para desarrollar la

i siguiente relación entre mediciones y predicciones del modelo:

j y = - 0 . 8 5 9 + 1.032x

i donde y = predicciones del modelo y x = mediciones. Conclu imos que había una buena

S concordancia entre los dos debido a que el intercepto fue aproximadamente igual a 0 y la

i pendiente, a 1. Recalcule la regresión pero ahora use la aproximación matricial para

i estimar los errores estándar de los parámet ros. Después emplee esos errore s para desa-

j rrollar los intervalos de confianza y úselos para hacer declaraciones probabilísticas con

i respecto a la bondad del ajuste.

¡] Solución. Los datos se pueden escribir en un formato matricial para regresión lineal

j simple como:

[Z] =

Se puede entonces usar la transposición y multiplicación de la matriz para generar las

ecuaciones normal como

"1 10 " 8.953

1 16.3 16.405

1 23 22.607

\Y) =

1 50 49 .988

[ÍZ]T[Z]] {A} = {ÍZ]T{Y}}

15 548.3 f 552.741

548.3 2219 1.21 _ [ 22421.43

La inversión de la matriz puede ser usada para obtener la pendiente e intercepto como

[A}= [íZflZ]]

0 .688414 -0 .01701

-0 .01701 0 .000465

{[Z]r{Y)

552.741

22421.43

- 0 . 8 5 8 7 2

1.031592



I

ii

Ii,i

1 7 . 4 FORMA GENERAL LINEAL

ruK/ v u N i m t »

wmmwww»De esta manera, el in tercepto y la pendie nte se determinan c om o a Q = — 0.85872y <7 »1031592respectivamente. Estos valores a su vez se pueden usar para calculnr el er rorestándar del es t imado como s /x =0863403Este valor puede u ti l izarse junto con loselementos d iagonal de la matr iz inversa para calcular los er rores es tándar de los coef icientes ,

La esta díst ica , í̂ n-i nec esa ria pa ra un 95% de int erv alo de co nfi an za co n n — 2

= 15 — 2= 13 g rado s de liberta d pue de ser dete rmi nada de una tabla estad ística omed iante sof tware. Usam os una función de Exce l , TINV, para obte ner el valor adec uado ,como en

=TNV00513)la cual da un valor de 2160368 Las ecuac iones (1729)y (1730) se pueden en toncesusar para calcular los in tervalos de conf ianz a com o

a0

= -085872 ±2160368(0716372= -085872 ±1547627=[-240634 0688912]

a , = 1031592±2160368(0018625)

=1031592±0040237=[09913551071828]Obse rve que los valores deseado s (0 para el in tercepto y pendiente, y 1 para el

in tercepto) es tán dentro de los in tervalos . Sobre la base de este anális is podr ía mo s hacerlas s iguientes declaracione s con respecto a la pendien te: tenemo s fuertes fundam entospara pensar que la pendiente de la l ínea de regresión real es tá dentro del in tervalo de0991355 a 1071828 Debido a que 1 está dentro de este in tervalo , tenemos tambiénbases fuer tes para creer que el resultado sopor ta la concordancia e ntre las me dicion es yel modelo . Co mo ce ro está dentro del in tervalo del in tercepto , una decla ración s imilar sepuede hacer con respecto del in tercepto .

Lo anter ior es una in troducción l imitada al enr iqueced or tema de la inferencia es tad ís t ica y su relación con la regresión . Hay much as más que están lejos del alcance de estel ib ro . Nue stra pr incipal in tención ha s ido mostrar el poder del enfoque m atr icial al ajustegeneral l ineal por mínim os cuadrado s. Usted deber ía consultar algunos de los excele nteslibros sobre el tem a (por ejem plo, Dr ape r y Smith 1981 para informac ión ad icio nal.Ade má s, deber ía observar qué paquetes de sof tware y l ibrer ías pueden gene rar ajustes deregresión por mínimos cuadrados junto con información relevante para la es tadís t icainferencial . Exploraremos algunas de esas capacidades cuando se descr iban esos paquetes al final del cap ítul o 19



4 9 6 REGRESIÓN PO R MÍNIMOS CUADRADOS

1 7 . 5 R E G R E S I Ó N N O L I N E A L

Hay muchos casos en ingenier ía donde modelos no l ineales deben ser ajustados condatos . En el contexto actual , esos modelos se def inen como aquellos que t ienen dependencia no l ineal de sus parámetros . Por ejemplo ,

Esta ecuación no puede ser manejada de acuerdo con el formato general de la ecuación

Com o se h izo con los mínim os cuadra dos l ineales , la regresión no l ineal se basa enla determ inación de los valores de los parámetros q ue minim izan la suma de los cua dra

dos de los res iduos. Sin embargo, para el caso no l ineal , la so lución debe proc eder en unaforma iterativa.

E l método de Gauss-Newton es un algor i tmo para minimizar la suma de los cuadrados de los res iduos en tre datos y ecuacione s no l ineales . El concepto clave que resalta latécnica es que una expansión por ser ie de Taylor se usa para expresar la ecuación nolineal or ig inal en una forma l ineal aproximada. Entonces, la teor ía de mínimos cuadrados se puede usar para obtener nuevas est imaciones de los parám etros que se mueve n enla d irección de minimizar el res iduo.

Para i lustrar cómo se hace esto , pr imero se puede expresar de manera general larelación entre la ecuación no l ineal y los datos c omo

y¡ - f {x¡; ao, a u . . . , am) + e¡

donde y ¡ = valor medido de la var iable dependiente, /"^, ; a 0, ax,..., a m) = ecuac ión que

es una función de la var iable independie nte x¡ y una función no l ineal de los parámetrosa 0, at,..., am, y e, = er ror aleator io . Por conveniencia, es te mod elo se pued e expresar deforma abreviada al omitir los parámetros ,

El mode lo no l ineal puede ser expandido dentro de una ser ie de Taylor alrededor devalores de paráme tro y reducido después de las pr imeras der ivadas. Por ejemplo , para uncaso de dos parámetros

(17.31)

(17.23) .

(17 .32)

/(*;),-+, - f(x¡)j +df(x¡\

da0

ja 0 +

a / t a ) .(17 .33)

donde j = son los valores iniciales, j + 1 = predicción, Aa 0 = « 0 j+i ~~ a

oj> y =aij+i ~ aij - ^ e esta forma, hemos l inear izado el modelo or ig inal con respecto a losparám etros . La ecuació n (17 .33) se pued e susti tu ir en la (17 .32) para obtener

o en forma matr icial [compárela con la ecuación (17 .24)] ,

\i>\ = [ / , J | A / \ ¡ , ( / : ¡ (17 .34)



donde [Z¡\ es la matriz de las derivadas parciales de la función evaluada en el valor ini

cial y,

[Z, ] =

dfi/dao dfi/dcn3 / 2 / 3 a o df2/dai

df„/da 0 dfn/daxdonde n = número de datos y df¡/dak = der ivada parcial de la función con respecto al k-

ésimo parámetro evaluado en el / -és imo punto . El vector {£>} contiene las d iferenciasentre las med iciones y los valores de la función,

y\ - f(x\)

>'2 - f(x 2)

¡D} =

y„ - f(x„)

y el vector {AA} contiene los cambios en los valores de los parámetros ,

ArtoAa\

{AA} =Aa„i se aplica la teor ía de mínimo s cuadra dos l ineales a la ecuación ( 17 .34) resulta en lass iguientes ecuaciones normal [ recuerde la ecuación ( 17 . 25) ] :

[[Z J]T[Z J}]{AA} = {[Z /]

T{D}}17 .35)

Así, el procedimiento consis te en resolver la ecuación ( 17 .35) para {A/4}, la cual sepuede emplear para calcular valores mejorados para los parámetros , como en

a0

j+i = O Q J + AaQ

Este proced imie nto se repite hasta que la so lución converge (es decir , hasta)

\S„h =k,j+\100% ( 17 .36)

está por debajo de un cr i ter io de paro aceptable.



498 R E G R E S I Ó N P O R M Í N I M O S C U A D R A D O S

E J E M P L O 1 7 . 9 M é t o d o d e G a u s s - N e w t o n

I Enunciado del problema. Ajuste la función f(x;a0, «,) = a

0( 1 ) con los dalos:

X 0 ,25 0 .75 1.25 1.75 2.25

y 0.28 0 .57 0 .68 0 .74 0 . 79

Use los valores in iciales de a0

= 1.0 y a x =1.0 para los parámetros . Observe que paraestos valores la suma in icial de los cuadrad os de los res iduos es 0 .024 8.

Solución. Las der ivadas parciales de la función con respecto a los parámetros son

da o = 1

= aQje

(E17.9 .1)

(E17.9 .2)

Las ecuaciones (E17.9 .1) y (E17.9 .2) se pueden usar para evaluar la matr iz

[Z0] =

02212 0194705276 0354307135 0358108262 0304108946 02371

Esta matr iz mult ip l icada por su traspuesta resulta en

'23193 09489"09489 04404Z][Z]

mi de nueve

[[Z][Zo]p=

la cual de nuevo se puede inver t ir para dar

36397 -78421-78421 191678

El vector {D } consis te en las d iferencias entre las mediciones ylas predicciones delm o d e l o ,

0280212 ' 0058805705276 0042406807135 — , -0033507408262 -0086207908946 -01046

{D } =

Ésta es mult ip l icada por [Z0]T para obtener

¡Z ()]T{D }

o í-01533-00365



17 . 5 REGRES IÓN N O L IN EAL mEl vector {AA} es entonces calculado al resolver la ecuación (17 .35) para

A A =- 0 . 2 7 1 4

0 . 5 0 1 9

la cual puede ser agregada al parám etro in icial supuesto para dar

I«o

I1. 0

1.0 + I- 0 . 2 7 1 4

0 . 5 0 1 9 i 0 . 7 2 8 6

1.5019 1

Así, las es t imaciones mejoradas de los parámetros son a 0 = 0.7286 y a x = 1.5019. L o snuevos parám etros resultan en una sum a de los cuadrad os de los res iduo s igual a 0 .0242.La ecu ación (17 .36) se puede usar para calcula r £ Q y e x igua l a 37 y 33 % , respec t ivamente . El cálculo podr ía repetirse hasta que esos valores es tén abajo del cr i ter io de paropreescr i to . El resultado f inal es a 0 = 0.79186 y a, = 1 .6751. Estos coef icientes dan unasuma de los cuadrado s de los res iduos de 0 .000662.

Un problema potencial con el método de Gauss-Newton como se ha desarro lladohasta ahora es que las der ivadas parciales de la función pue den se r d if íci les de evaluar .En consecuencia, muchos programas de computadora usan d iferentes ecuaciones paraaproximar las der ivadas parciales . Un método es

dfi ^ f(x¡; a 0,...,a k+ Sa k, . . . , a m ) - f{x¡\a 0, ..., a k,..., am)

da k 5a k

donde 8 = per turbación f raccional pequeña.El método de Gauss-Newton t iene una var iedad de o tros posib les defectos:

1 . Pue de converger con lenti tud .2 . Puede oscilar ampliamente; es decir , cambia en forma continua de d irección.3. Puede no converger .

Se han desarro llado modif icacion es del méto do (Booth y Peterson, 195 8; Har t ley , 1961)para remediar los defectos .

Adem ás , aunque hay var io s p roced imien tos exp resamen te d iseñados para r eg res ión ,uno más general es usar ru tinas de optimización no l ineal como las descr i tas en la par tecuatro . Para hacer es to , se hace una suposic ión de los parám etros , y se calcula la suma de

los cuadrados de los res iduos. Por ejemplo , para la ecuación (17 .31) es to se podr ía calcula r como

Entonces, los parámetros se podr ían ajustar de manera s is temática para minimizar S r

mediante técnicas de búsqueda descr i tas previamente en el capítu lo 14 . I lustraremos elmo do para hacer es to cuando desc r ibamo s las aplicaciones del sof tware al final del capí-

(17 .38)

tu lo 19 .



R E G R E S IÓ N P O R M Í N I M O S C U A D R A D O S

P R O B L E M A S

17.1 Dados los datos

0.90 1.42 1.30 1.55 1.631.32 1.35 1.47 1.95 1.661.96 1.47 1.92 1.35 1.051.85 1.74 1.65 1.78 1.712.29 1.82 2.06 2.14 1.27

determine a) la media, b) la desviación estándar, c) la varianza,d) el coeficiente de variación y e) el intervalo de confianza al95% para la media.

17.2 Construya un histograma para los datos del problema 17.1.Use un rango de 0.6 a 2.4 con intervalos de 0.2.17.3 Con los datos

15 6 18 21 26 28 3239 22 28 24 27 27 33

2 12 17 34 29 31 3845 36 41 37 43 38 26

determine a) la media, b) la desviación estándar, c) la varianza,d) el coeficiente de variación y é) el intervalo de confianza al90% para la me dia./) Construya un histograma. Use un rango de0 a 55 con incrementos de 5. g) Suponga que la distribución esnormal y que su estimación de la desviación estándar es válida,

calcule el rango (esto es , los valores inferior y superior) que abarquen el 6 8% de las lecturas. Determine si esto es una estimaciónválida para los datos en este problema.17.4 Use la regresión por mínimos cuadrados para ajustar a unalínea recta a

X l 3 7 10 12 13 16 18 20

y 4 5 6 5 8 7 6 9 12 11

Junto con la pendiente y el intercepto, calcule el error estándardel estimado y el coeficiente de correlación. Grafique los datosy la línea de regresión. Después repita el problema, pero ahorahaga la regresión de x contra y (es decir, cambie las variables).

Interprete sus resultados.17.5 Use regresión por m ínimos cuadrados para ajustar una línea recta a

X 5 6 10 14 16 20 2 2 28 28 36 38

30 22 28 14 22 16 8 8 14 0 4

Junto con la pendiente y el intercepto, calcule el error estándardel estimado y el coeficiente de correlación. Grafique los datosy la línea de regresión. Si alguien hizo una medición adicionalde x = 5, y = 5, ¿podría usted esperar, con base en un asegura

miento visual y en el error estándar, que la medición fue válid a ofallida? Justifique su conclusión.17.6 Use regresión por mínimos cuadrados para ajustar una línea recta a

x 2 3 4 7 8 9 5 5

y I 9 6 5 K) 9 ñ 2 3

a) Junto con la pendiente y el intercepto, calcule el errorestándar de la estimación y el coeficiente de correlación.Grafique los datos y la línea recta. Valide el ajuste.

b) Recalcule a), pero ahora use regresión polinomial para ajustar los datos a una parábola. Compare los resultados conaquellos de a) .

17.7 Ajuste un modelo de razón de crecimiento saturado con

X 0.75 2 2.5 4 6 8 8.5

y 0.8 1.3 1.2 1.6 1.7 l . í3 1.7

Grafique los datos y la ecuación. Encuentre el error estándar.17.8 Ajuste los datos del problema 17.7 a una ecuación de potencias. Grafique los datos y la ecuación y encuentre el errorestándar.17.9 Adecúe los datos del problema 17.7 a una parábola. Grafiquelos datos y la ecuación y encuentre el error estándar.

17.10 Ajuste los siguientes datos a una ecuación de potenciasX 2.5 3.5 5 6 7.5 10 12.5 15 17.5 20

/ 7 5.5 3.9 3.6 3.1 2.8 2.6 2.4 2.3 2.3

Grafique y contra x junto con la ecuación de potencias.17.11 Adecué los siguientes datos a un modelo exponencial

X 0.4 0.8 1.2 1.6 2.0 2.3

y 7 5 0 1 0 0 0 1 4 0 0 2 0 0 0 2 7 0 0 3 7 5 0

Grafique los datos y la ecuación en papel estándar así como engráfico logarítmico. Analice sus resultados.17.12 Ajuste los datos del problema 17.11 a una parábola.

Grafique los datos y la ecuación.17.13 Con los datos

X 5 10 15 2 0 25 30 35 40 45 50

y 16 25 32 33 38 36 39 40 4 2 4 2

use regresión por mínimos cuadrados para ajustar a) a una línearecta, b) a una ecuación de potencias, c) a una ecuación de razónde crecimiento saturado y d) a una parábola. Grafique los datosjunto con todas las curvas. ¿Es superior alguna de las curvas? Sies así, justifique.



PROBLEMAS 5 0 1

17.14 Use regresión lineal múltiple para ajustar

x, 0 1 1 2 2 3 3 4 4

0

y 15 18 12.8 25 .7 20 .6 35.0 29.8 45.5 40.3

Calcule los coeficientes, el error estándar del estimado y el coeficiente de correlación.17.15 Use regresión lineal múltiple para ajustar

0 0 1 2 0 1 2 2 1

x 2 0 2 2 4 4 6 ó 2 1

y 15 19 12 11 24 22 15 5 19

Calcule los coeficientes, el error estándar del estimado y el coeficiente de correlación.

17.16 Use regresión no lineal para ajustar los siguientes datos auna parábola

y

0 .07 5 0.5 1.2 1.7 2 .0 2 .3

6 0 0 8 0 0 1 2 0 0 1 4 0 0 2 0 5 0 2 6 5 0 3 7 5 0

17.17 Use regresión no lineal para ajustar los datos del problem a17.13 a una ecuación de razón de crecimiento saturado.17.18 Recalcule el ajuste por regresión de los problemas a) 17.4y b) 17.13 mediante el procedimiento de la matriz. Estime loserrores estándar y desarrolle intervalos de confianza al 90% parala pendiente y el intercepto.17.19 Desarrolle, depure y pruebe un subprograma en ya sea unlenguaje de alto nivel o en lenguaje macro de su elección paraimplementar regresión lineal. Entre otras cosas: a) A gregue co

mentarios para documentar el código y b) determine el errorestándar y el coeficiente de determinación.17.20 Use el software de Métodos Numéricos TOOLKIT pararesolver los problemas a) 17 .4, b) 17.5,c) 17.6, d) 17.9 ye) 17.12.



C A P Í T U L O 1 8

Interpolación

Usted a menudo habrá tenido la opor tunidad de est imar valores in termedios entre datosprecisos . El método más común que se usa para es te propósito es la in terpolación del

polinomio. Recuerde que la fórmula general para un polinomio de n-ésimo orden es

f(x) = a 0 + aix + a2x2 + •••+ a nx" (18.1)

Para n + 1 puntos , hay uno y só lo ún polinom io de orden n que pasa a través de todos lospuntos . Por ejemplo , hay sólo una l ínea recta (es decir , un polinomio de pr imer orden)que conecta dos puntos (vea la figura 18 .1a) . De man era s imilar , única men te una pa rábola conecta un conjunto de tres puntos (ver figura 18.1Z>). Interpolación polinomial

consis te en determinar el único polinomio de n-ésimo orden que ajuste n + 1 pun tos .Este polinomio entonces proporciona una fórmula para calcular valores in termedios.

Aunque hay uno y só lo un polinomio de n-ésimo orden que ajusta n + 1 pun tos ,exis te una var iedad de formatos matemáticos en los cuales es te polinomio puede expresarse. En este capítu lo descr ib iremos dos al ternativas que son muy adecuadas para laimp lemen tac ión en compu tado ra : lo s po l inomios de Newton y de La g range .

F I G U R A 1 8 . 1

Ejemplos de interpolación polinomial: a] de primer orden [lineal) conectando dos puntos,b) de segundo orden (cuadrática o parabólica) enlazando tres puntos y c) de tercer orden(cúbica) conectando cuatro puntos.

a) b) c)



18.1 DIFERENCIA DIVIDIDA DE N E W TO N PARA LA INTERP OLAC IÓN 5 0 3

1 8 . 1 D I F E R E N C I A D I V I D I D A D E N E W T O N P A R A L AI N T E R P O L A C I Ó N D E P O L I N O M I O S

Como se mencionó antes , ' ex is te una var iedad de formas al ternativas para expresar unain terpolación de un polinomio. La diferencia dividida de Newton para la interpolación

de polinomios está entre los mo delos m ás populares y ú ti les . Ante s de presentar la ec uación general , in troduciremos la pr imera y segunda vers ión debido a su in terpretaciónvisual s imple.

1 8 . 1 . 1 I n t e r p o l a c i ó n l i n e a l

El modo más s imple de in terpolación es conectar dos puntos con una l ínea recta. Estatécnica, l lamada interpolación lineal, se ilustra en forma gráfica en la figura 18.2. Mediante tr iángulos semejantes ,

/ l ( - r ) ~ / ( S o ) = / ( * ! ) - / ( * » ) s

X — XQ X\ — Xola cual se pued e reordenar para dar

/ , ( * ) = / ( % > + / ( j C l ) ~ / ( X o ) (* - %> (18-2)

que es una fórmula de interpolación lineal. La no tac ión f(x) des igna que es una

in terpolación de polinomios de pr imer orden. Observe que además de representar la

F I G U R A 1 8 . 2Ilustración gráfica de interpolación lineal. Las áreas sombreadas indican los triángulossemejantes usados para obtener la fórmula de interpolación lineal [véase la ecuación! 1 8.2)]

AQ X X-\ X



rrvrcKruLACTurN

pendiente de la l ínea que conecta los puntos , el término \f(xx) — f(x

0)]/(x

x— x 0 ) es una

aproximac ión por d iferencia d iv id ida f in i ta de la pr imera der ivada [ recuerde la ecuac ión(4.17)] . En general , cuanto más pequeño sea el in tervalo de datos , mejor será la aproximación. Esto se debe al hecho de que, en tanto el in tervalo d isminuya, una funcióncontinua se aproximará mejor por una l ínea recta. Esta caracter ís t ica se demuestra en els iguiente ejemplo .

EJ EM PLO 18 . 1 I n terp ola c ión l inea l

Enunciado del problema. Estime el logar i tmo natural de 2 mediante in terpolación

lineal . Pr ime ro , realice el cálculo por in terpolación entre ln 1 = 0 y ln 6 = 1 .791759.Después, repita el procedimiento , pero use un in tervalo más pequeño de ln 1 a ln 4(1 .386294) . Obse rve que el valor real de ln 2 es 0 .6931472.

Solución. Usaremos la ecuación (18 .2) y una in terpolación l ineal para ln(2) de x0 = 1

a X] = 6 para dar

1 791759 - 0/ , ( 2 ) = 0 + — (2 - 1) = 0 . 3 5 8 3 5 1 9

6 — 1

que representa un er ror de e t = 4 8 . 3 % . Con el in tervalo más pequeño de, x0 = 1 a x

x=

4 se obtiene

/ , ( 2 ) = 0 +

L 3 8

4

6 ^ - ° (2 - 1 ) = 0 . 4 6 2 0 9 8 1

Así, con el in tervalo más cor to se reduce el er ror relat ivo porcentual a e t — 3 3 . 3 % .Ambas in terpolaciones se muestran en la f igura 18.3, junto con la función real .

F I G U R A 1 8 . 3Dos interpolaciones lineales para estimar ln 2 . Obse rve cómo el intervalo más pequeñoproporciona una mejor estimación.

f (x) = ln j f .

Es t imac iones l i nea les

l i l i l í



18.1 DIFERENCIA DIVIDIDA DE N E W TO N PARA LA INTERP OLA CIÓN 5 0 5

1 8 . 1 . 2 I n t e r p o l a c ió n c u a d r á ti c a

El er ror en el ejemplo 18.1 resulta de nuestra aproximación a una curva con una l ínearecta. Por consiguiente, una estrategia para mejorar la es t imación es in troducir algunacu rva tu ra en la l ínea que conec ta los pun tos . S i t i e s^ u iü o^ de ja s da to s es tán d ispon ib les ,es to puede r ea l iza r se con un po l inomio de segundo o rden ( tamb ién conoc ido comopolinomio cuadrático o parábola). Una forma en par t icular conveniente para es te propósito es

fi(x) = bo +bi(x- x 0) + b 2(x - x 0){x - x{) (18 .3)

Observe que aunque la ecuación (18 .3) parecer ía d ifer ir del polinomio general [véaseecuación (18 .1)] , las dos ecuaciones son equivalentes . Esto puede demostrarse al mult ip licar los términos de la ecuación (18 .3) para dar

f 2(x ) = b0

+ b\x - b]X 0 + b2x2 + b

2x

Qx\ - b

2xx

0- b

2xx \

o, agrupando términos,

f2(x ) = ao + a\x + a 2x

2

donde

a o = bo - btxo + h 2X{)X\

a\ = b\ — b2xu — b 2x\

a2

= b2

Así, las ecuac iones (18 .1) y (18 .3) son formulaciones al ternativas equivalentes del únic opolinomio de segundo orden que une los tres puntos .

Un procedim iento s imple puede usarse para determinar los valores de los coef icient e s . Para Z > 0 , la ecuación (18 .3) conx — X Q puede ser usada para calcular

b 0 = / ( * o ) (18 .4)

La ecua ción (18 .4) puede susti tu irse en la (18 .3) , la cual puede evaluarse en x = x x para

= f^)-f(xo)

X \ - XQ

Por ú lt imo, las ecuac iones (18 .4) y (18 .5) se pued en susti tu ir en la (18 .3) , la cual p ued e

evaluarse en x = x 2 y resolver (después de algunas manipulaciones algebraicas) para

/ ( * 2 ) - / ( * i ) / ( * i ) - / ( * o )

, xo -x\ xi - x0 (18 .6)

f> 2 = =

x2

- XoNote que, como fue el caso con la in terpolación l ineal , b x todavía representa la

pendiente de la l ínea que une los puntos x 0 y x x. Así, los pr imeros dos términos dela ecuación (18 .3) son equivalentes a la in terpolación l ineal de x 0 ax , , como se espec i f i có antes en la ecuación (18 .2) . El ú l t imo térm ino, b 2(x — x 0 ) (x — X j ) , in troduce la curvatura de segundo orden en la fórmula.



5 0 6 I N TERPOLACIÓN

E J E M P L O 1 8 . 2

ij

|¡í

Antes de i lustrar cómo usar la ecuación (18 .3) , debemos examinar la forma delcoef iciente b 2. Es muy s imilar a la aproximación por d iferencias d iv id idas f in i tas de lasegunda der ivada qu e se in trodujo antes en la ecuación (4 .24) . Así , la ecua ción (18 .3)comienza a manifestar una estructura que es muy s imilar a la ser ie de expansión deTay lo r . Es ta observac ión se r á ob je to de más exp lo rac ión cuando r e lac ionemos lo spolinom ios de in terpolación de New ton con la ser ie de Taylor en la sección 18.1 .4 . Peropr imero , desarro llaremos un ejemplo que muestra cómo se usa la ecuación (18 .3) parain terpolar entre tres puntos .

Interpolación cuadrát ica

Enunciado de l problema. Ajuste los tres puntos usados en el ejemplo 18.1 a unpolinomio de segundo orden:

xo = 1 f(xo) - 0

Xl=4 / ( x j ) = 1.386294

x2 = 6 / ( J C 2 ) = 1.791759

Use el polinomio para evaluar ln 2 .

Solución. Aplicando la ecuación (18 .4) se obtiene

b0 = 0

La ecuación (18 .5) da

F I G U R A 1 8 . 4Uso de interpolación cuadrática para estimar ln 2 . Para comparación se incluye también lainterpolación lineal de x = 1 a 4 .

2

0

1

x5



18.1 DIFERENCIA DIVIDIDA DE NE WTON PARA LA INTERPOLA CIÓN 507

y con la ecuación (18.6) se obtiene1.791759 - 1.386294

0.4620981

b2 = 6 - 4 ^ = -0.05187 316—1

Sustituyendo estos valores en la ecuac ión (18.3) se obtiene la fórmula cuadrática

f2(x) = 0 + 0.462098l(x - 1) - 0.0518731 (x - ])(x - 4)

que puede evaluarse en x — 2 para

/ 2 (2) = 0.5658444

la cual representa un error relativo de et = 18.4%. Así, la curvatura introducida por la

fórmula cuadráti ca (véase la figura 18.4) mejora la interpolación comparada con el re

sultado que se obtiene mediante líneas rect as en el ejemplo 18.1 y figura 18.3.

18 .1 .3 For ma general de la interpolaciónde polinomios de Newton

El análisis anterior puede ser generalizado para ajustar un polinomio de «-ésimo orden a

n + 1 datos. El polinomio de «-ésimo orden es

f„(x) = b0 + bi(x - x0) + • • • + b„(x - x0)(x -*,)•••(*- x„-i) (18.7)

Como se hizo antes con las interpolaciones lineales y cuadráticas, los puntos de los datos

evaluaban los coeficientes b0, b¡,..., b

n. Para un polinomio de «-ésimo orden se requiere

n + 1 puntos: [x 0 , / (x 0 ) ] , [x^fix^],. .., [xn,f(xn)]. Usamos estos datos y las siguientes

ecuaciones para evaluar los coeficientes:

bQ = f(x0) (18.8)

b\=f[xux0] (18.9)

bi = f[x z,xi, x„) (18.10)

bn = f[x„,x„-i, . . . , I l , A '

0

] (18.11)

donde las evaluaciones de la función puestas entre paréntesis son diferencias divididas

finitas. Por ejemplo, la primera diferencia dividida finita se representa por lo general

como

rr T / ( • * / ' ) ~ f (X

j ) / 1 0 n \

f[xj,Xj] = (18.12)X¡ - Xj

La segunda diferencia dividida finita, la cual representa la diferencia de las dos primeras

diferencias divididas, se expresa por lo general como



508 I NTERPOLACI ÓN

; * /

0 xo f ( x 0 ) -

ix\ f ( x , ) -

2 x 2 f ( x 2 | -

3 * 3 f[x 3h

P r i m a r a

: f [ x , , x 0 ] -

: f [ x 2 , x , ]

: f [ x 3 . x 2 ] "

S e g u n d o

| . x„]

Ta reero

| X | , X , , X ,„

F I G U R A 1 8 . 5

Ilustración gráfica de la naturaleza recursiva de las diferencias divididas finitas.

, r , f[Xi,Xj] - f[.Xj,X K]

f \x ¡ ,x ,,x ¡Á = - -

J 1 x¡-x k (18 .13)

En forma similar, la n-ésima diferencia dividida finita esj .

r, fí xni xn — \ > • • • i x \ \ — / [ x „ _ i , X „ _ 2 , . . . , XQ ] (W l\A \

f[x„,X„-u...,XUXo] = - - - - - ( . 1 8 . 1 4 )

X„ — XQ

Estas d iferencias pueden usarse para evaluar los coef icientes en las ecuaciones (18 .8)hasta la (18 .11) , las cuales entonces se susti tu irán en la ecuac ión (18 .7) para obtene r elpolinomio de in terpolación

/«(*) = ñ xo) + (x- * 0 y i * i » * o ] + ( x~ x0)(x - x^x^ xu x0]

+ -+ (x - Xo)(x - x, ) ... (x - x n_ 1Y[x n,x n_ l,...,x 0] (18 .15)

que es conoc ido como polinomio de interpolación por diferencias divididas de Newton.

De be obse rvarse que no es necesar io que los datos u ti lizados en la ecuación (18 .15) seanigualmente espaciados o que los valores de la abscisa deban estar en orden ascendente,com o se i lustra en el siguiente ejemplo . También, observe c óm o las ecuaciones (18 .12) a(18.14) son recursivas (es decir , las d iferencias de órdene s mayo res se calculan al tom ardiferencias de orden menor, véase la figura 18.5). Esta propieda d será aprovechada cuand odesarro lle mos u n programa ef iciente en la computadora dentro de la sección 18.1 .5 paraimp lemen ta r e l método .

EJ EM PL O 18 . 3 I n terp o la ndo p o l i nomi os med ia nte la d i fe rencia d i v i d i d a d e Newt o n

Enunciado del problema. En el ejemplo 18.2, los datos e n x 0 = \,x x = 4 yx 2 — 6 seuti l izaron para es t ima r ln 2 con una parábola . Ahora, ag regando un cu ar to punto (x 3 = 5;f(x 3) = 1 .609438] , calcule el ln 2 con una in terpolación del polinomi o de New ton detercer orden.

Solución. El polinom io de tercer orden u ti l izando la ecuación (18 .7) con n — 3, es

/ l ( . v ) = b{) + l> l(X - X T ) ) + l>2(.\ - .«o) (.V - .V|) | /0,(.\ .loK.v • -.»| )(\ .v>)



1 5.1 UirCKCrNUA U l V I V I U r t u c P I S Y T I W I ^

F IGURA 18 .6Uso de la interpolación cúbica para estimar ln 2.

Las primeras diferencias divididas para el problema son [véase la ecuación (18.12)]

1.386294-0f[X],X 0] =

f[x2,Xi] =

f[X3,X2] =

4 - 0

1.791759 - 1.386294

6 - 4

1.609438 - 1.386294

5 - 6

= 0.4620981

0.2027326

0.1823216

Las segundas diferencias divididas son [véase la ecuación (18.13)]

0.2027326-0.4620981f[X2,Xi,X 0] =

f[XT,,X 2,Xi

6 - 1

0.1823216-0.2027326

5 - 4

= -0.05187311

= -0.02041100

La tercera diferencia dividida es [véase la ecuación (18.14) con n = 3]

- 0 . 0 2 0 4 1 1 0 0 - (-0.05187311)

5 - 1= 0.007865529

Los resultados par a/[ *,, x0],f[x2, X ^ X Q ] yf[x3, x2, xx, x0] representan los coeficientes 6,,

b 2 y ¿>3 de la ecuación (18.7). Junto con b0 =f(x0) = 0.0, la ecuación (18.7) es



310 INTERPOLACIÓN

fr(x) = 0 + 0.4620981 ( J I - 1) - 0.05187311 (x I )(.v - 4)

+ 0.007865529(x - l)(x - 4)( x - 6)

que puede usarse para ev aluar/ ^) = 0.6287686, el cual representa un error relativo de

e¡ = 9.3%. El polinomio cúbico completo se muestra en la figura 18.6.

18 .1.4 Errores al interpolar polinomios de Newton

Observe que la estructura de la ecuación (18.15) es similar a la serie de expansión de

Taylor en el sentido de que se agrega términos en forma secuencial para capturar el

comportamiento de alto orden de la función en turno. Estos términos son diferencias

divididas finitas y, así, representan aproximaciones de las derivadas de orden mayor. Por

consiguiente, como ocurrió con la serie de Taylor, si la verdadera función subyacente es

un polinomio de w-ésimo orden, el polinomio sujeto a interpolación de n-ésimo con base

en n + 1 puntos dará resultados exactos.

También, como fue el caso con la serie de Taylor, puede obtenerse una formulación

para el error de truncamiento. Recuerde de la ecuación (4.6) que el error de truncamien

to para la serie de Taylor podría expresarse por lo general como

f (n+l)(t\

(n+ 1)1(4.6)

donde ¿j está en alguna parte en el intervalo x¡ a x, + 1 . Para una interpolación de n-ésimo

orden, una relación análoga para el error es

R„ r"+1(£)(n + 1)!

(x - x 0 ) ( x - x t ) • • • (x -x„) (18.16)

donde £ está en alguna parte en el intervalo que contiene la incógnita y los datos. Para

esta fórmula que habrá de usarse, la función en turno debe ser conocida y < iferenciable.

Por lo común éste no es el caso. Por fortuna, una formulación alternativa ei tá disponible

y no requiere conocimiento previo de la función. Más bien, usa una diferí ncia dividida

finita para aproximar la derivada (n + l)-ésima,

R„ = fíx,x„,x„^i , XQ](X - X0)(X - X\) • • • (X - X„)

dondef[x, x n,x nA,..., x0)] es la (n + 1 )-ésima diferencia dividida finita. D :bido

ecuación (18.17) contiene la incógnita/(x), no puede resolverse para el err >r.

go , si se dispone de un dato adicional f(x n+1), la ecuación (18.17) pue<eestimar el error, como en

R„ = f[x„+u x „ , . . . , x0](x - x0)(x - x,) • • • (x - x„) (18.18)

E J E M P L O 18.4 Es t imación del er ror pa ra el pol inomio de Newton

Enunciado del prob lema. Use la ecuación (18.18) para estimar el error para la

¡ interpolación del polinomio de segundo orden del ejemplo 18.2. Use los datos adiciona-

! les/(x 3 ) = /(5) = 1.609438 para obtener sus resultados.

(18.17)

a que la

Sin embar-

usarse para



• • » .* n u i n r w i v i v n v e r w H m w m i u s ' n i m i w n O f f ¥ 1 7

F I G U R A 1 8 . 1 0Una ilustración visual del racional detrás del polinomio de Lagrange. Esta figura muestra uncaso de segundo orden. Cada uno de los tres términos en la ecuación ( 18 .2 3 ) pasa através de uno de los puntos y es cero en los otros dos. La sumatoria de los tres términos, portanto, es un polinomio único de segundo orden f 2|x) que pasa exactamente a través de lospuntos.

Observe que, com o en el caso del método de New ton, la vers ión de Lagrange t ieneun error estimado de [véase ecuación (18.17)]

n

R n = f[x,x„,x„-i,... ,x0] Y\(x -x¡)

De este modo, si se t iene un punto adicional en x — xn+l

, se puede obtener un errorestimado. Sin embargo , como no se emplea n las diferencias divididas f initas co mo partedel algoritmo de Lagrange, esto se hace muy ocasionalmente.

Las ecuac iones ( 18 . 20) y ( 18 . 21) s e programan de manera muy s imple paraimplementarse en una computadora. La figura 18.11 muestra el pseudocódigo que sepuede emplear para este propósito.

En resumen, para casos donde el orden del polinomio es desconocido, el método deNewton tiene ventajas en el conocimiento que proporciona respecto al comportamientode las fórmulas de diferente orden. Ad em ás, el error estimado e xpuesto por la ec uació n(18,18) se integra de manera usual en el cálculo del polinomio de Newton debido a que



5 1 8 I NTERPOLACI ÓN

FUNCTION Lagrng(x, y, n, x)

sum = 0

DO i = 0 ,n

product - y

DO = 0 ,n

IF i # j THEN

product = product*(x - Xj ) / (x ; - Xj)

ENDIF

EWDO

sum = sum + product

EWDO

Lagmg = sum

END Laqrftq

F I G U R A 1 8 . 1 1Pseudocódigo para ¡mplementar la interpolación de Lagrange. Este algoritmo seacondiciona para calcular una sola predicción de n-ésimo orden, donde n + 1 es el númerode datos.

la es t imación em plea una d iferencia f in i ta (véase el ejemplo 18.5) . De esta man era, paracálculos explorator ios , a menudo se pref iere el método de Newton.

Cuando se ejecuta só lo una in terpolación, las formulaciones de Lagrange y Newtonrequieren un notable esfuerzo co mputa cional . Sin emba rgo, la vers ión de Lagra nge es unpoco más fácil de programar . Debido a que no requiere de cálculos y almacenaje dediferencias d iv id idas , la forma de Lagrange se usa a menudo cuando el orden del polinomiose conoce a priorí.

EJ EM PL O 1 8 . 7 I n terp o la c ión d e La g ra ng e usa nd o la comp ut a d ora

Enunciado del problema. Podem os usar el algor i tmo de la f igura 18 .11 para e s tudiarun anális is de tendencia para un problema que se relaciona con nuestio conocido casodel paracaid is ta en caída. Suponga que desarro llamos cier ta instrumentación para medula velocidad del paracaid is ta . Lo s datos de med ición obtenid os para un caso de prueba enpar t icular son

T i e m p o ,s

V e l o c i d a d m e d i d av , e m / s

1 8003 2 3105 3 0 9 07 3 940

13 4 755

í Nue stro problema es es t imar la velocidad del paracaid is ta en t — 10 s para tener lasmediciones fal tantes entre t — lyt= 13 s . Estam os conscientes que el com por tam ientode la in te rpo lac ión de po l inomios puede se r inesperado . Po r tan to , cons t ru i r emospolinomios de 4 , 3, 2 y 1 órdenes para comparar los resultados.



18 . 3 COEF I C I ENTES D E UN POL I NOMI O D E INTERPOLACI ÓN

| F I G U R A 1 8 . 1 2] Gráfica s que exponen a) cuarto orden, b) tercer orden, c) segundo orden¡ y d) interpolaciones de primer orden.

Solución . El algor i tm o de Lagra nge se usa para constru ir polino mio s de in terpolaciónde cuar to , tercer , segundo y pr imer orden.

El polinomio de cuar to orden y los datos de entrada se pueden graf icar como semue stra en la f igura 18 .12a. Es evidente a par t ir de esta gráf ica que el valor es t imad o dey e n x = 10 es mayor que la tendencia g lobal de los datos .

Desd e la f igura 18 .12b hasta la 18.12d se muestran las gráf icas de los resultados d elos cálculos para las in terpolaciones de los polinomios de tercer , segundo y pr imer or den. Se observa que mientras más bajo sea el orden, menor será el valor es t imado de lavelocidad en t = 10 s . Las gráf icas de la in terpolación de polinomios indican que lospolino mio s de al to orden t ienden a sobrepasar la tendencia de los datos . Esto sugiere q uelas vers iones de pr imer o segundo orden son las más adecuadas para es te anális is detendencia en par t icular . Debe recordarse, que debido a que tratamos con datos incier tos ,la regresión d e hecho podr ía ser la má s adecuada.

El ejem plo anter ior i lustró que los polinomios de al to orden t ienden a ser ma l condicionado s; es to es , t ienden a ser al tam ente sensib les a los er rores de redondeo. El m ism o

problema se aplica a la regresión con polinomios de orden super ior . La ar i tmética dedoble precis ión ayuda algunas veces a d isminuir el problema. Sin embargo, en tanto elorden aum ente , habrá un punto para el cual el er ror de redondeo in ter fer irá con la hab il idad para in terpolar mediante el procedimiento s imple que se abordó en este punto .

1 8 . 3 C O E F I C I E N T E S D E U N P O L I N O M I O P E I N T E R P O L A C I Ó N

Aunque el polinomio de Newton y el de Lagrange son adecuados para determinar valores in termedios entre puntos , no proporcionan un polinomio conveniente de la formaconvenc iona l




f(x) =a 0 + a xx + a 2x2 + •••+ a nx

n (18 .24)

Un método d irecto para calcular los coef icientes de este polinomio se basa en elhecho de que n + 1 puntos se requieren para determinar los n + 1 coef icien tes . Así , sepuede usar ecuaciones algebraicas l ineales s imultáneas para calcular las a. Por ejemplo ,suponga que usted desea calcular los coef icientes de la parábola

f(x) = a0

+ a\x + a2x 2 (18 .25)

Se requiere tres punto s: [x 0,f(x 0)] , [ X ] , / ( x j ) ] y [x 2,f(x 2)] .Cada uno se pued e susti tu ir enla ecuación (18 .25) para dar

f(x0) — a

0+ aix

0+ a

2x\

f(x\) = a 0 + a\x\ + a 2x\ (18 .26)

f(x 2) = a 0 + a\x 2 + a 2x\

De esta manera, las x son los puntos conocidos y las a las incógnitas . Como hay elmismo número de ecuaciones que incógnitas , la ecuación (18 .26) se podr ía resolver conun método de el iminación de la par te tres .

Debe observarse que el procedimiento anter ior no es el método d isponible más ef iciente para determinar los coef icientes de una in terpolac ión de polinom ios. Press y cois .(1986) proporcionan un anális is y códigos de cómputo para procedimientos más ef icaces . Cualquiera que sea la técnica empleada, se debe tener precaución con el orden. Lossis temas como el de la ecuación (18 .26) es tán notor iamente mal condicionados. Ya seanresueltos con un método de el iminación o con un algor i tmo ef iciente, los coef icientesresultantes puede n ser al tamente inexactos , en par t icular para n grandes . C uando se u sanpara una in terpolación subsecuente a menudo se obtiene resultados er róneos.

En resumen, s i usted se in teresa en determinar un punto in termedio , emplee lain terpolación de Newton o Lagrange. Si lo que desea es determinar una ecuación con laforma de la (18 .24) , l imítese a polinomios de orden infer ior y ver if ique sus resultadoscu idadosamen te .

1 8 . 4 I N T E R P O L A C I Ó N I N V E R S A

Como su nomenclatura lo implica, los valores de f(x) y x en la mayo r ía de los contextosde in terpolación son las var iables dependientes e independientes , en forma respectiva.En consecuencia, los valores de las x están t íp icamente espaciados en forma uniforme.

Un ejemplo s imple es una tabla de valores tabulados para la función f(x) = l/x,

X1 2 3 4 5 6 7

1 0.5 0 . 3333 0.25 0 .2 0 . 1 6 6 7 0 . 1 4 2 9

Ahora suponga que usted debe usar los mismos datos , pero que se le ha dado unvalor para f(x) y debe determinar el valor correspondiente de x. Por ejemplo , para losdatos anter iores , supong a que se le p ide dete rmina r el valor de x que co r r esponda a f(x)— 0.3 . Para es te caso , com o la función está d isponible y se puede man ejar en forma fácil ,la respuesta correcta se determina directamente com o A ; = 1/0.3 = 3.3333.



punt os en una secuencia d iferente. La figura 18.9 mue stra los resu ltad os para el cas o deinvertir el orden de los datos originales; esto es, x 0 = 3 .5 , x { = 2.5 , x 3 = 1.5, y asisucesivamente. C om o los puntos in iciales para es te caso se hallan ce rcanos y e spacia dosen cualquier lado de ln 2 , el er ror d isminuye mucho más rápido que para la s i tuaciónor ig inal . En el segundo término, el er ror se redujo a menos de e, = 2% . Se pod r ía emplear o tras combinaciones para obtener d iferentes velocidades de convergencia.

El ejemplo anter ior i lustra la impor tancia de la selección de los puntos base. Comose podr ía in tu ir en forma obvia, los punto s deber ían estar centrado s alrededo r y tan c ercacom o sea posib le de las incógnitas . Esta observación es tam bién sop or tada por un anális is d irecto de la ecuación para es t im ar el er ror [vea ecuación (18 .17)] . Si supo nem os que

la d iferencia d iv id ida f in i ta no var ía mucho a lo largo del rango de datos , el er ror esproporcional al producto: (x — x 0) (x — x¡).. .(x — x n). Por lógica, los puntos base máscercanos son a x, la magnitud más pequeña de este producto .

1 8 . 2 I N T E R P O L A C I Ó N D E P O L I N O M I O S P E L A G R A N G E

L a interpolación de polinomios de Lagrange es s imp lemen te una r e fo rmu lac ión de lpolinomio de Newton que evita el cálculo por d iferencias d iv id idas . Se puede expresarde manera conc isa com o

/»(*) = 2 £,(*)ft*i> (18.20)

i= 0

donde

i W =n^ / (18.2Dy=o X Í - X J

donde I I designa el "producto de". Por ejemplo , la vers ión l ineal (n = 1) es

fi(x) = ?—^f(xo) + í—^-f(xx) (18.22)x 0 - xi x, - x 0

y la vers ión de segundo orden es

' ( i - x i ) ( x - x 2 ) - ( x - x 0 ) ( x - X2

)

f2(x) = - — r / ( x 0 ) + - — r / ( * i )( x 0 - xO(x o - x 2 ) {x\ - x 0)(xi - x 2)

, (x-xo)(x-Xl)m fi X2 ) ( 1 8 . 2 3 )

(x 2 - x 0 ) ( x 2 - X i ) "

La ecuación (18 .20) se puede der ivar de manera d irecta a par t ir del polinomio deNewton (véase cuadro 18 .1) . Sin embargo, el racional resaltado de la formulación deLagrange se puede captar de manera d irecta al darse cuenta que cada término L ¡(x) será




1 en x = x, y 0 en todos los otros pun tos de la mue stra (véa se la figura 18.10). De estamanera, cada producto L¡(x)f(x¡) toma el valor de f(x¡) en el punto de muestra x¡ . E nconsecuencia, la sumator ia de todos los productos designados para la ecuación (18 .20)es el único polinom io de n-ésimo orden qu e pasa de manera exa cta a través de todos losn + 1 puntos .

EJ EM PLO 1 8 . 6 I n terp ola c ión d e p o l i nomi os d e La g ra ng e

Enunciado del problema. Use una in terpolación del polinomio de Lagrange de pr ime r y segundo orden para evaluar ln 2 con base en los datos dados en el ejem plo 18.2 :

xo = 1 f(xo) = 0

x , = 4 / ( x 0 = 1.386294

x2=6 f(x

2) = 1.791760

Solución. El polinomio de pr imer orden [véase ecuación (18 .22)] se puede usar paraobtener la es t ima ción en x = 2 ,

i / , ( 2 ) = Y~^°+ y 1 - 3 8 6 2 9 4 = 0 . 4 6 2 0 9 8 1

De manera s imilar , el polinomio de segundo orden se desarro lla como [véase ecuación1 (18 .23)]

( 2 - 4 ) ( 2 - 6 ) ( 2 - l ) ( 2 - 6 )/ 2 ( 2 ) = - 0 + - 1 .3 8 6 2 9 4

( l - 4 ) ( l - 6 ) ( 4 - l ) ( 4 - 6 )( 2 - l ) ( 2 - 4 )

+ - - 1 .791 760 = 0 .5658444( 6 - l ) ( 6 - 4 )

Como se esperaba, ambos resultados concuerdan con los que se obtuvieron antes al usarla in terpolación para el polinomio de Newton.

Cuadro 18 .1 Obtención de la forma de Lagrange a partir directamentede la interpolación del polinomio de Newton

El polinomio de interpolación de Lagrange se puede obtener demanera directa a partir de la formulación del polinomio deNewton. Haremos esto únicamente para el caso del polinomiode primer orden [véase ecuación (18.2)]. Para obtener la formade L agrange, reformulamos las diferencias divididas. Por ejemplo, la primera diferencia dividida,

fl i /(xi) - /(x0)/[x,,xo] = (B18.1.1)xi - xose puede reformular com o

la cual es referida com o orma simétrica. Al sustituir la ecuación(B18.1.2) en (18.2) se obtiene

/i(x) = /(x0) + — — f t x i ) + ——/(x0)x\ - x0 x0 - X i

Por último, al agrupar términos similares y simplificar se tienela forma del polinomio de Lagrange,

X — X l x — Xof x ( x ) = -/(x0) + V(xi)xo - X l X\- X Q

(B18 .1 .2 )



18.5 COMENTARIOS ADICIONALES 521

Tal problema se conoce como interpolación inversa. Para un caso más complicado,

usted debe intentar cambiar los valoresf(x) y x [es decir, sólo grafique x contra f(x)\ y

use un procedimiento como la interpolación de Lagrange para determinar el resultado.

Por desgracia, cuando usted invierte las variables no hay garantía de que los valores

junto con la nueva abscisa [los/(x)] sean espaciados uniformemente. De hecho, en mu

chos casos, los valores serán "agrandados" . Es decir, tendrán la apariencia de una escala

logarítmica con algunos puntos adyacentes muy agrupados y otros muy dispersos. Por

ejemplo, para f(x) = 1/x el resultado es

0.1429 0.1667 0.2 0.25 0.33 33 0.5 1

X 7 6 5 4 3 2 1

Tal espaciamiento no uniforme sobre la abscisa a menudo lleva a oscilaciones en el

resultado de la interpolación de polinomios. Esto puede ocurrir aun para polinomios de

orden inferior.

Una estrategia alterna es ajustar un polinomio de n-ésimo orden, fn{x), a los datos

originales [es decir, con f(x) contra x]. En la mayoría de los casos, como las x están

espaciadas de manera uniforme, este polinomio no estará mal condicionado. La respues

ta a su problema entonces toma en cuenta la determinación del valor de x que haga que

este polinomio sea igual al dado por f(x). Así, ¡el problema de interpolación se reduce a

un problema de raíces!

Por ejemplo, para el problema antes descrito, un simple procedimiento podría ser

ajustar los tres puntos a un polinomio cuadrático: (2, 0.5), (3, 0.3333) y (4, 0.25). El

resultado sería

f2(x) = 1.08333 - 0.375.* + 0.041667*

2

La respuesta al problema de interpolación inversa para determinar lax correspondiente a

f(x) = 0.3 sería equivalente a la determinación raices de

0.3 = 1.08333 - 0.375A- + 0 .041667A : 2

Para este caso simple, la fórmula cuadrática se puede usar para calcular

_ 0.375 ± v/C-0.375) 2 - 4(0.041667)0.78333 _ 5.704158X ~ 2(0.041667) ~ 3.295842

Así, la segunda raíz, 3.296, es una buena aproximación del valor real de 3.333. Si se

desea una exactitud adicional, se podría emplear un polinomio de tercer o cuarto orden

junto con uno de los métodos para la localización de raíces de la parte dos.

18 .5 C O M E N T A R I O S A D I C I O N A L E S

Antes de proceder con la siguiente sección, se debe mencionar dos temas adicionales:

interpolación con datos igualmente espaciados y extrapolación.

Como, ambos polinomios, el de Newton y Lagrange, son compatibles con datos

espaciados en forma arbitraria, usted se preguntaría por qué hemos puesto el caso es

pecial de datos igualmente espaciados (véase cuadro 18.2). Antes de la llegada de las

computadoras digitales, esas técnicas tenían gran utilidad para interpolación a partir de



5 2 2 INTERPOLACIÓN

C u a d r o 1 8 . 2 I N T E R P O L A C I Ó N C O N D A T O S I G U A L M E N T E E S P A C I A D O S

Si los datos están igualmente espaciados y en orden ascendente,entonces la variable independiente tiene los valores de

x\ = x 0 + h

X2 = AO + 2h

donde lo que resta es lo mismo que en la ecuación (18.16). Estaecuación es conocida como fórmula de Newton o fórmula hacia

adelante de Newton-Gregory. Se puede simplificar más al definir una nueva cantidad, a:

x„ = xo + nh

donde h es el intervalo, o tamaño de paso, entre los datos. Conesta base , las diferencias dividid as finitas se pueden expresar enforma conc isa. Por ejemplo, la segunda diferencia dividida haciaadelante es

/[A'O, A|, X 2] =

f(Xl) - f(Xj) _ / ( A l ) ~ / ( X 0 )

A 2 — X l X| — XQ

X2 -X0la cual se puede expresar como

/ ( x 2 ) - 2 / ( x , ) + / ( x 0 )[ X 0 , Al , A 2 ]2h 2

(B 18.2.1)

ya que x x — x 0 = x 2 — x { = (x 2 — x 0)/2 = h . Ahora recuerde

que la segunda diferencia hacia adelante es igual al numeradorde la ecuación (4.24)

A 2 / ( x o ) = / ( A ' 2 ) - 2 / ( x l ) + / ( x 0 )

Por tanto, la ecuación (B 18.2.1) se puede representar como

A 2 / ( x 0 )/[ X 0,X|,X2 ] =

o, en general

/ [ X o , X | , . . . , x „ ]

2\h 2

= A"/(x 0 )

n\h "(B18.2.2)

Mediante la ecuación (B l 8.2.2), podemos expresar el polinomiode interpolación de Newton [véase ecuación (18.15)] para el casode datos igualmente espaciados como

Esta definición se puede usar para desarrollar las siguientes expresiones simplificadas pára los términos en la ecuación (B18.2.3):

x — A O = ahx — A Q — h = ah • •h = h(a-\)

x — Ao — (n — l)h = ah — (n — \)h = h(a — n + 1)

las cuales se sustituyen en la ecuación (Bl 8.2.3) para dar

A 2 fíx0)fn (x) = / ( x 0 ) + A/(x 0 )a + i) "a(a - 1)

donde

2 !, A"/(x 0 )

H -i ¡ —a ( a - 1) • • ( a - n + 1) + R „

(B 18.2.4)

R n = , , , \ , V + W - D (« - 2) • • • (a - «)[n + 1)!

Esta notación concisa tendrá utilidad en nuestra derivación y enel análisis de error de las fórmulas de integración en el capítulo21 .

Además de la fórmula hacia adelante, están disponibles lasfórmulas hacia atrás y central de Newton-Gregory. Para información adicional con respecto a interpolación para datos igualmente espaciados se puede consultar a Carnahan , Luther y Wilkes(1969).

r / x s, , A / ( X 0 )/ , ( X ) = / ( X Q ) + ; ( X - X 0 )

h

+A 2 / ( A 0 )

2\h 2(x - x 0 ) ( x - xo - h)

+ ---+^¡^(x-xQ)(x-x0-h)nlh"

•••\x-x0-(n-1)h) + R „ (B18.2.3)



I i o. i uircRCNClA DIVIDIDA DE NE W TO N "PARA LA INTERPO LACIÓN

IMSolución. Recuerde que en el ejemplo 18.2 la interpolación del polinomio de segundoj orden proporcionó una estimación de f2(2) = 0.5658444, que representa un error de

\ 0.6931472 - 0.5658444 = 0.1273028. Si se hubiera conocido el valor real, como es

" común que suceda, la ecuación (18.18), junto con el valor adicional en x3, pudo usarse

para estimar el error, como en

! Ri = / f e , x 2 , X\,XQ\(X - x0)(x - x\)(x - x2)

io

! R2 = 0.007865529(JC - l)(x - 4)(x - 6)

í donde el valor para la diferencia dividida finita de tercer orden es como la que se calculó

antes en el ejemplo 18.3. Esta relación puede evaluarse en * = 2 paraR

2 = 0.007865529(2 - 1)(2 - 4)(2 - 6) = 0.0629242

la cual es del mismo orden de magnitud que el error real.

Del ejemplo anterior y de la ecuación (18.18), debe estar claro que el error estimado

para el polinomio de n-ésimo orden es equivalente a la diferencia entre el (n + l)-ésimo

orden y la predicción de n-ésimo orden. Es decir,

Rn = fn +l(x) ~ fn(x) (18.19)

En otras palabras, el incremento que se agrega al caso del n-ésimo orden para crear el

caso de (n + l)-ésimo orden [es decir, la ecuación (18.18)] se interpreta como unaestimación del n-ésimo orden de error. Esto puede verse con claridad al reordenar la

ecuación (18.19) para dar

f„+l{x) = f„(x) + Rn

La validez de este procedimiento es afirmada por el hecho de que la serie es fuertemente

convergente. Para tal situación, la predicción del (n + l)-ésimo orden debería ser mucho

más cercana al valor real que la predicción del n-ésimo orden. En consecuencia, la ecua

ción (18.19) conforma nuestra definición estándar de error al representar la diferencia

entre la verdadera y una aproximación. Sin embargo, observe que mientras todos los

otros errores estimados para procedimientos iterativos introducidos hasta ahora se deter

minaron como una predicción actual menos una previa, la ecuación (18.19) representa

una predicción futura menos una presente. Esto significa que para una serie en convergencia rápida, el error estimado de la ecuación (18.19) podría ser menor que el verdade

ro . Esto representaría una calidad muy poco atractiva si el error estimado fuera a emplearse

como un criterio de paro. Sin embargo, como será expuesto en la siguiente sección, la

interpolación de polinomios de alto orden son muy sensibles a errores en los datos (es

decir, están mal condicionados). Cuando se emplean para interpolación, a menudo se

obtienen predicciones que divergen en forma significativa del valor verdadero. Si se

prevén tales errores, la ecuación (18.19) es más sensible para dicha divergencia. Como

tal, es más valioso para la clase de análisis de datos exploratorios para lo cual el polinomio

de Newton es el más adecuado.




1 8 . 1 . 5 A l g o r i t m o de c ó m p u t o p a r o la i n t e r p o l a c i ó nd e l p o l i n o m i o de N e w t o n

Tres propiedades hacen que la in terpolación del polinomio de Newton sea muy atractivapara las aplicaciones en la compu tado ra :

1. C o m o en la ecuación (18 .7) , se puede desarro llar de manera secuencial vers iones de

orden mayor con la adición de un so lo té rmino a la s igu ien te ec uac ión de ordeninferior. Esto facilita la evaluación de a lgunas de las vers iones de diferente orden en

e l m i s m o p r o g r a m a . Tal c a p a c i d a d es en espec ia l va l io sa cuando el o r d e n del

po l inomio no es conocido a priori. Al ag regar nuevos té rminos en forma secuencial ,podemos de te rminar cuándo se a lcanza un pun to de d isminuc ión de r eg reso (es

decir , cuando la adición de t é rminos de orden super ior ya no mejo ra de m a n e r asignif icativa la estimación, o en cier tas s i tuaciones de h e c h o la aleja) . Las ecuac io nes para es t imar el error que se ana l izan en el pun to 3 son úti les para v isualizar un

criterio objetivo para identificar este punto de t é rminos d isminu idos .

2 . Las diferencias divididas finitas que constituyen los coeficientes del polinomio [ecuaciones (18.8) hasta (18.11)] se pueden ca lcu la r de manera e f icaz . Es decir , como en

la ecuación (18 .14) y la figura 18.5, se usa diferencias del orden infer ior para calcular las de alto orden. Por me dio de la in fo rmación de te rminada an tes , los coef icientes se pueden ca lcu la r de manera ef iciente. El a lgo r i tmo en la f igura 18.7 con t ieneun esquema semejan te .

3 . El error es t imado [véase la ecuac ión (18 .18 ) ] puede ser muy s imp le de incorporaren un a lgo r i tmo de cómpu to deb ido a la manera secuenc ia l en la cual se construye la

pred icc ión .

F I G U R A 18 . 7

Un algoritmo para el polinomio de interpolación deNewton escrito en pseudocódigo.

5U3R0UTINE Newtlnt (x, y, n, xi, ymt, ea)

LOCAL fdd„¡n

DO i = 0,n'

fddu>=y¡

END DO

DO = 1,n

DO i = 0,n- j

fdd¿ = (fddM¡t i - fdd^x^-xj

ENDDOEND DO

xterm = 1

yint0 = fdd0¡0

DO order =1n

xterm = xterm * (x¡ - x o r e t e r _ , )

yintZ = yintolder.1 + fdd0Order * xterm

END order

END Newtlnt



18 .1 D IF ERENCIA D IV IDIDA DE N EW T O N P ARA L A INT ERP OL AC IÓN 91 3

Todas las caracter ís t icas anter iores pueden aprovecharse y ser incorporadas en unalgor i tmo general para impleme ntar el polinom io de New ton (véase la f igura 18 .7). Observe que el algor i tmo consis te en dos par tes : el pr imero determina los coef icientes opar t ir de la ecuación (18 .7) ; el segundo establece las predicc iones y su er ror asociado. Lauti l idad de este algor i tmo se demuestra en el ejemplo s iguiente.

EJ EM PL O 1 8. 5 Es t imación del erro r para determinar e l orden adecuado de in terpolación

¡ Enunciadod el problema. Después de incorporar el er ror [véase la ecuación (18 .18)] ,; u t i l ice el algor i tmo de cómputo que se mues tra en la f igura 18 .7 y la s iguiente informa -j ción para eva lua r /(x ) = Inxenx = 2 :

i x

f(x ) = ln x

0 14 1.3862944

ó 1.7917595

5 1.6094379

3 1.09861231.5 0 .40546412 .5 0 .9 1 6 29 0 73.5 1.2527630

Solución. Los resultados al emplear el algor i tmo de la f igura 18 .7 para obtener unasolución se muestran en la f igura 18 .8 . El er ror es t imado, junto con el er ror real (con

| F I G U R A 1 8 . 8

Lo resuados de un pogamacon base en e agomo de la fgua18.7,para evauar ln2.

NUMBER OF P01NTS? 8

X ( 0 ) , y ( 0 ) - ? 1 ,0

X ( 1 ) . y ( 1 ) - ? 4 , 1 3 8 6 2 9 4 4

X ( 2 ) , y ( 2 ) - ? 6 , 1 . 7 9 1 7 5 9 5

X ( 3 ) . y ( 3 ) - ? 5 , 1 . 6 0 9 4 3 7 9

X( 4 ) , y ( 4 ) - ? 3 , 1 . 0 9 8 6 1 2 3

X( 5 ) . y ( 5 ) - ? 1 .5 , 0 . 4 0 5 4 6 4 1 1

X( 6 ) , y ( 6 ) - ? 2 . 5 . 0 . 9 1 6 2 9 0 7 3

X( 7 ) , y ( 7 ) - ? 3 . 5 , 1 . 2 5 2 7 6 3 0

INTER POLA TION AT X - 2

ORDER F ( X ) ERROR

0 0 . 0 0 0 0 0 0 0 . 4 6 2 0 9 8

1 0 . 4 6 2 0 9 8 0 . 1 0 3 7 4 6

2 0 . 5 6 5 8 4 4 0 . 0 6 2 9 2 4

3 0 . 6 2 8 7 6 9 0 . 0 4 6 9 5 3

4 0 . 6 7 5 7 2 2 0 . 0 2 1 7 9 2

5 0 . 6 9 7 5 1 4 - 0 . 0 0 3 6 1 6

6 0 . 6 9 3 8 9 8 - 0 . 0 0 0 4 5 9

7 0 . 6 9 3 4 3 9



514I NTERPOLACI ÓN

base en el hecho de q ue ln 2 = 0 .6931472) , se i lustran en la f igura 18 .9 . Obs erve q ue eler ror es t imado y el real son s imilares y que su concordanc ia mejo ra en tanto se a ume nteel orden. A par t ir de estos resultados se puede concluir que la vers ión de quinto orden dauna buena est im ación y que los término s de orden super ior no resaltan en forma s ignif icativa la predicción.

I Este ejercicio tam bién i lustra la impo r tancia de colocar el orden de los puntos . Por¡ ejemplo , hasta la es t imac ión del tercer orden, la razón de mejora e s lenta debido a que\ los puntos que se agregaron (en x = 4, 6 y 5) están distantes y a un lado del punto de! análisis en x = 2. La est im ación de cuar to orden muestra algo de mejora ya que el nuevo\ punto en x — 3 está cercano a la incógnita . Sin embargo, la d isminución más dramática\ en el er ror es tá asociada con la inclusión del térmi no de quinto orden median te los datos

\ en x = 1 .5 . No sólo está es te punto cercano a la incógnita , s ino que tambié n se halla en¡ el lado opuesto de la mayor ía de los o tros puntos . En consecue ncia, el er ror se reduce ai casi un orden de mag nitud .j El s ignif icado de la posición y secuen cia de los datos se pue de tam bién i lustrar alj usar los mism os datos para obtener una est imac ión para el ln 2 , pero conside rando los

i F I G U R A 1 8 . 9j Porcentaje de errores relativos para la predicción de ln 2 como una función! del orden de la interpolación polinomial.

Error i



F I G U R A 1 8 . 1 3Ilustración de la divergencia posible de una predicción extrapolada. La extrapolación sebasa en el ajuste de una parábola a través de los primeros tres puntos conocidos.

tab las con argumentos igualmente espaciados. De hecho, una estructura computacionalconocida como tabla de d iferencias d iv id idas fue desarro llada para facil i tar la imple-men tación de esa s técnica s . (La f igura 18 .5 es un ejemplo de esa tabla. )

S in embargo , como las fó rmu las son subcon jun tos de esquemas de Newton yLagrange compatib les con computadora y debido a las muchas funciones tabulares d isponibles , com o subrutinas de l ibrerías , la necesidad para las vers iones igua lme nte e spaciadas ha d ism inuido. A pesar de esto , por su relevancia las hem os inclu ido en este te maen las ú lt ima s par tes de este l ibro . En par t icular , son necesar ias pa ra obtener fórm ulas dein tegración numérica que emplean de manera t íp ica datos igualmente espaciados (véaseel capítu lo 21) . Com o las fórmulas de in tegración num érica t ienen relevancia en la so lución de ecuaciones d iferenciales ordinar ias , el mater ial del cuadro 18 .2 t iene tambiénsignif icado en la par te s iete.

Extrapolación es el proceso de est im ar un valor de f(x) que se t iene fuera del rangode los puntos base conocidos, x 0, J C , , . . ., x n (véase la figura 18.13). En una secciónanter ior , mencionamos que la in terpolación más exacta es usualmente obtenida cuandolas incógnitas es tán cerca del centro de los puntos base. Obviam ente, es to no se cum plecuando las incógnitas se encuentran fuera del rango y , en consecuencia, el er ror enextrapolació n pued e ser mu y grande. Co mo se i lustra en la figura 18.13, la naturaleza deextremo s abier tos de la extrapolación representa un paso en la incógnita , ya que el proceso extiende la curva más al lá de la región conocida . Com o tal , la curva real podr ía confacil idad d ivergir de la predicción. Por tanto , se debe tener extremo cuidado al realizarejercicios donde surja un caso que se deba extrapolar.



INTERPOLACÓNI N T E R P O L A C I Ó N S E G M E N T A R I A

-^n la sección anterior, se usó polinomios de w-ésimo orden para interpolar entre n +1

¿Jatos. Por ejemplo, para ocho puntos se puede derivar un perfecto polinomio de séptimo¿>rden. Esta curva podr ía capturar todas las curvaturas (al meno s hasta e incluso la séptica derivada) sugeridas por los punto s. Sin em bargo , hay casos en los que estas funcione spu ed en llevar a resultados er róneos d ebido a er rores de redonde o y puntos lejan os. U n

sfGURA 1814(Jna representación visual de una situación donde las segmentarias son interpolacionespolinomialesde orden superior. La función que habrá de ajustarse pasa por un incrementojubito en x = 0. Los incisos a) a c) indican que el cambio abrupto induce oscilaciones al¡nterpolar polinomiales. En contraste, como las curvas se limitan a tercer orden contransiciones suaves, la segmentaria cúbica d ) proporciona una aproximación mucho más¿iceptable.

a)

f(x),

Xb)

f(x) iJ> 9 9 9

— , »-

J •0 xd)



18 . 6 I NTERPOLACI ÓN SEGMENTARI A 525

procedimiento al ternativo es aplicar polinomios de orden infer ior a subconjuntos de datos . Tales polinomios conectores son l lamados funciones segmentarias.

Por ejemplo , curvas de tercer orden empleadas para conectar cada par de datos sonl lamadas segmentarias cúbicas. Esas funciones se pueden constru ir de tal forma que lasconexiones entre las ecuaciones cúbicas adyacentes resultan v isualmente suaves. Sobrela superf icie, podr ía parecer que la aproximación de tercer orden de las segmentar iassería inferior a la expresión d e séptimo orden . UstecJ se pregun taría por qu é una segm entariapodr ía ser s iempre prefer ib le.

La f igura 18 .14 i lustra una s i tuación donde una segme ntar ia se com por ta mejo r q ueuna polinomial de orden super ior . Este es el caso donde una función es por lo generalsuave pero conlleva un camb io abrupto e n algún lugar a lo largo de la región de in terés .

¡ í El increm ento de paso expuesto en la f igura 18 .14 es un ejemplo extrem o de tal cam bio

y s irve para i lustrar es te punto .De la figura 18 .14a hasta la 18 .14c se ilustra cómo un polinom io de orden su per iort iende a formar una curva a través de oscilaciones brusca s en la vecindad co n un cam biosúbito . En contraste, la segmentar ia también conecta los puntos , pero debido a sus cambios l imitados de tercer orden, las oscilaciones se mantienen a un mínimo. Como tal , lasegmen ta r ia u sua lmen te p ropo rc iona una ap rox imación super io r de l compor tamien to delas funciones que t ienen cambios locales y abruptos .

El concepto d e la segmentar ia se or ig ina de la técnica de d ibujo con una cin ta de lgada y flexible (llamada curvígrafo) para d ibujar curva s suaves a través de un conjunto depuntos . El proceso se expone en la f igura 18 .15 para una ser ie de cinco pasadores (datos ) . En esta técnica, el d ibujante coloca un papel sobre una mes a de madera y golpea losclavos o pasadores en el papel (y la me sa) e n la ubicació n de los datos . Una curva suaveresulta al en trelazar la cin ta entre los pasa dores . De aquí que se haya adoptad o el nom bre

de "segmentar ia cúbica" para los polinomios de este t ipo .En esta sección, se usarán pr imero funciones l ineales s imples para in troducir algu

nos conceptos básicos y problemas asociados con la in terpolación segmentar ia . Entonces obtendremos un algor i tmo para el ajuste de datos con segmentar ias cuadráticas . Porúlt imo, presentamos mater ial sobre la segmentar ia cúbica, la cual es la vers ión más comú n y ú ti l en la práctica de la ingenie r ía .

F I G U R A 1 8 . 1 5\i La técnica de dibujo.al usar

una segmentaria para

dibujar curvas suaves a

través de una serie de

. puntos. Obse rve cóm o en

los puntos extremo, la

segmentaria trata de

', ende rezars e. Esto es

J cono cido como uña

1 segmentaria "natural" .




1 8 . 6 . 1 S e g m e n t a r i a s l i n e a l e s

La conexión má s s imple entre dos puntos es por medio de una l ínea recta. Las segmentar iasde pr imer orden para un grupo de datos ordenados pueden def in irse como un conjuntode funciones l ineales ,

f(x) = f(xo) +m 0(x - x 0 )

/ ( x ) = / ( x i ) +n n(x - x i )

x o < x '< xi

X i £ x < Xi

f(x) = / ( x „ - j ) + m„-\(x - x „ _ i ) x „ _ i < x < x„

donde m¡ es la pendiente de la l ínea recta que conecta los puntos :

/(*/ + !) - f(x¡)

x¡+i — x¡(18 .27)

Estas ecuac iones se pueden usar para evaluar la función en cualq uier punto entre x 0

y x n para localizar pr imero el in tervalo dentro del cual es tá el punto . Después se usa laecuación ad ecuada pa ra determinar el valor de la función dentro del in tervalo . El mé todoes obviame nte idéntico al de la in terpolación l ineal .

EJ EM PL O 1 8 . 8 Seg ment a r i a s d e p r i mer o rd en

Enunciado de l problema. Ajuste los datos de la tabla 18 .1 con segmentar ias de pr imer orden. Evalúe la función enx = 5 .

Solución. Se puede usar los datos para determinar las pendientes entre puntos . Porejemplo , para el in tervalo x = 4 .5 a x = 7 la pendiente se puede calcular mediante laecuación (18 .27) :

2.5 - 1

7 - 4 . 5= 0.60

Las pendie ntes para los o tros in tervalos se pued en calcular y las segmenta r ias res ultantes

de prim er orden se grafic an en la figura 18.16a. El valor en x = 5 es 1.3.

T A B L A 1 8 . 1Datos para ser ajustadoscon funciones segmentarias.

3.04. 57.09. 0

f ( x )

2.51.02.50.5



Una inspección visual a la figura 18.16a indica que la principal desventaja de las

segm entar ias de pr imer orden es que no son suaves. En esenc ia, en los puntos d onde dossegmentar ias se encuentran ( l lamado nodo), la pendiente cambia de forma abrupta. Entérm inos formales , la pr imera der ivada de la función es d iscontinua en esos puntos . Estadef iciencia se resuelve al usar segmentar ias polinomiales de orden super ior que aseguren suavidad en los nodos al igualar der ivadas en esos puntos , como se analiza en las iguiente sección.

1 8 . 6 . 2 S e g m e n t a r i a s c u a d rá t ic a s

Para asegurar que las der ivadas m-ésimas son continuas en los nodos, se debe usar unasegmentar ia de al menos m + 1 orden. A menudo se usan con más f recuencia en lapráctica los polinomios de tercer orden o segmentar ias cúbicas para asegurar der ivadas

F I G U R A 1 8 . 1 6Ajuste segmentario de un conjunto de cuatro puntos, a) Segmentaria lineal, b) segmentariacuadrática y c] segmentaria cúbica, se gráfica también con una interpolación polinomialcúbica.

Segmentaria de

primer orden

10 x

Segmentaria de

segundo orden

0 I 1 1 | L J L

b)

' o

Segmentaría

cúbica

\Interpolación

cúbica

f




continuas de primer y segundo orden. Aunque las derivadas de tercer orden y mayor©podrían ser discontinuas cuando se usa segmentarias cúbicas, usualmente no puededetectarse en forma visual y en consecuencia son ignoradas.

De bido a que la obtención de segm entarias cúbicas está algo involucrada, la hem*escogido en una sección subsecuente. Hemos decidido primero ilustrar el cor.cípto e jinterpolación segmentaria mediante polinomios de segundo orden. Esas "segmentáis jcuadráticas" tienen primeras derivadas co ntinuas en los n odo s. Aunqu e las segmentarascuadráticas no aseguran segund as derivadas igu ales en los nodo s, «'"•en muy b¡en pía |demostrar el procedimiento general en el desarrollo de segmentarias de o. jen ¡superior.

El objetivo de las segm entarias cuadráticas es derivar un polino mio de segundo orden para cada intervalo entre datos. El polinomio para cada intervalo se puede represen

tar de manera general co mo

f , ( x ) = c n x2

+ b i x + c i ( 1 8 . 2 8 )

La figura 18.17 ha sido incluida para ayudar a clarificar la notación. Para n + 1 datos(i = 0, 1, 2 n) existen n intervalos y, en consecuencia, 3n constantes desconocidas(las a,byc) por evaluar. Por tanto, se requieren 3« ecuacio nes o co ndicio nes para evaluar las incógnitas. Éstas son:

1 . Los valores de la función de polinomios adyacentes deben ser iguales a los nodosinteriores. Esta condición se puede representar como

a¡ -xx}_

x+ + c,-i = / ( * , _ ! ) ( 1 8 -2 9 )

a¡xf_ + b¡x¡-[ + c¡ = /(* , •_ i ) (18 .30)

para / = 2 a n. Com o sólo se usa nodos interiores , las ecuaciones (18 .29) y (18 .30)proporcionan cada una n — 1 cond iciones del total de 2 n — 2 .

F I G U R A 1 8 . 1 7

Notación usada paraderivar segmentariascuadráticas. Observe quehay n intervalos y n + 1datos. El ejemplo mostradoes para n = 3 .

a,)? + bfX+ c-

Intervalo 1•« •

Intervalo 2* •

Intervalo 3

_ j 1 , 1 1* 0 * 1 * 2 * 3 *

; = 0 / = 1 1-2 / = 3



1 8 . 6 * I N T E R P O L A C I Ó N S E G M E N T A R I A 1 9 9

2 . La s primera y última funciones deben pasar a través de lo s puntos extremo. Eagrega dos ecuaciones adicionales:

a\xl + b ix 0 + c\ = f(x0) (18.31)

a„ x 2n + bnx„ +c„- f(xn) (18.32)

para un total de 2n — 2 + 2 = 2n condiciones.

3 . Las primeras derivadas en lo s nodos interiores deben ser iguales. La primera deri

vada de la ecuación 18.28 es

f\x) = 2ax + b

Por tanto, la condición se puede representar de manera general como

2a¡-\x,-\ + b¡-[ = 2 a¡Xi-\ + b ¡ (18.33)

para i = 2 an. Esto proporciona otras n — 1 condiciones para un total de 2n + n —

1 = 3« — 1. Como se tiene 3n incógnitas, se tiene una condición corta. A menos

que tengamos alguna información adicional con respecto a las funciones o sus deri

vadas, debemos tomar una selección arbitraria para calcular de manera exitosa las

constantes. Aunque hay un número de elecciones diferentes que se pueden tomar,

seleccionamos la siguiente:

4 . Suponga qu e en el primer pu nto la segunda derivada es cero. Como la segunda

derivada de la ecuación 18.28 es 2a¡, esta condición se puede expresar matemática

mente como

«i = 0 (18.34)

La interpretación visual de esta condición es que los dos primeros puntos se conec

tarán con una línea recta.

E J E M P L O 18.9 Segmentarias cuadráticas

j Enunciado del problema. Ajustar por medio de segmentarias cuadráticas los mismos

! datos que se usaron en el ejemplo 18.8 (véase tabla 18.1). Use los resultados para calcu-

\ lar el valor enx = 5.

| Solución. Para este problema, se tienen 4 datos conn— 3 intervalos. Por tanto, 3(3) =

\ 9 incógnitas por ser determinadas. Las ecuaciones (18.29) y (18.30) dan 2(3) — 2 = 4i condiciones:

1 2 0 . 2 5 a i + 4 . 5 & 1 + c i = 1.0

I 20 .25(3 2 + 4.5Í72 + c 2 = 1.0

j 49a 2 + 7 / 7 2 + c 2 = 2.5

J 49 a 3 + 7fo3 + c 3 = 2.5

Evaluando a las funciones primera y última en los valores inicial y final se agregan 2

ecuaciones mes [véase ecuación (18.31)]:



5 3 0 INTERPOLACIÓN

9a! +3b> +c i = 2.5

y [véase ecuación (18.32)]

81a3 + 9¿> 3 +c 3 = 0.5

La continuidad de las derivadas crea un adicional de 3 — 1 = 2 [véase ecuación (18.33)]:

9Ú [ + b\ = 9a2 4- b2

14a 2 + ¿>2 = 14a3 + ¿>3

Por último, la ecuación (18.34) especifica que ax

— 0. Como esta ecuación especifica a¡de manera exacta, el problema se reduce a la resolución de ocho ecuaciones simultáneas.

Estas condiciones se pueden expresar en forma matricial como

4.5 1 0 0 0 0 0 0" 1

0 0 20.25 4.5 1 0 0 0 C\ 1

0 0 49 7 1 0 0 0 a 2 2.5

0 0 0 0 0 49 7 1 b22.5

3 1 0 0 0 0 0 0 C 22.5

0 0 0 0 0 81 9 1 a 3 0.5

1 0 - 9 - 1 0 0 0 0 0

0 0 14 1 0 - 1 4 - 1 0 C '3 0

Esas ecuaciones se pueden resolver mediante las técnicas de la parte tres, con los resultados:

« i = 0 b\ = — 1 c\ = 5.5

a2 = 0.64 b2 = - 6 . 7 6 c2 = 18.46

a 3 = - 1 . 6 ¿ 3 = 2 4 . 6 £-3 = - 9 1 . 3

las cuales pueden sustituirse en las ecuaciones cuadráticas originales para desarrollar la

siguiente relación para cada intervalo:

/, (.v ) = —jr + 5.5 3 . 0 < x < 4 . 5

f2(x) = 0 .64 A - 2 - 6.76J T + 18.46 4.5 < x < 7.0

/ 3(.v) = -1 . 6 x2

+ 24.6x - 91.3 7.0 < x < 9.0

Cuando se usa^, la predicción para x = 5 es, por tanto,

/ 2 (5) = 0.64(5)2 - 6.76(5) + 18.46 = 0.66

El ajuste total por segmentarias se ilustra en la figura 18.166. Observe que hay dos

desventajas que se alejan del ajuste: 1) la línea recta que conecta los primeros dos puntos

y 2) la segmentaria para el último intervalo parece oscilar demasiado. Las segmentarias

cúbicas de la siguiente sección no exhiben estas desventajas y, en consecuencia, son

i«<*molación seamentaria.



1 8 . 6 . 3 S e g m e n t a r i a » c ú b ic a s

El objetivo en las segmentarias cúbicas es obtener un pol inomio de tercer orden partícada intervalo entre los nodos, como en

f¡ ( J C ) = a/.v3 + b¡ x 2

+ c¡x + d-, ( 1 8 . 3 5 )

Así, para n + 1 datos (i — 0, 1, 2, . .. , n), existen n intervalos y, por consiguiente, 4«incógnitas constantes para evaluar. Com o co n las segmentarias cuadráticas, se requieren4 H condicio nes para evaluar las incógnitas. Éstas son:

1 . Los valores de la función deben ser iguales en los nodos interiores (2« — 2 condiciones) .

2 . La primera y última fun ciones de ben pasar a través de los puntos extremo (2 condiciones) .

3. Las primeras derivadas en los no dos interiores deben ser iguales (« — 1 condiciones) .4 . Las segundas derivadas en los nodo s interiores deben ser iguales (« — 1 condiciones) .5. Las segundas derivadas en los nodos extremo son cero (2 condiciones) .

La interpretación visual de la condición 5 es que la función se vue lve una línea recta en

los nodos extremo. La especif icación de tal condición extremo nos lleva a lo que se

denomina como segmentaria "natural". Se le da este nombre debido a que el dibujosegmentario se comporta en esta forma (véase f igura 18.15).

Si el valor de la segunda derivada en los nodo s extremo no es cero (es decir, existe algunacurvatura), esta información se puede usar de manera alterna para suministrar las dos

condiciones f inales .Los cinco t ipos de condiciones anteriores proporcionan un total de 4n ecuaciones

requeridas para resolver los An coeficientes. Mientras sea ciertamente posible desarro

llar segmentarias cú bicas en esta forma, presentamos una técnica alterna que requiere la

solución de sólo n — 1 ecuaciones . Aunque la obtención de este método (véase cuadro18.3) es algo menos directo que para las segmentarias cuadráticas, la ganancia en eficiencia bien vale el esfuerzo.

La derivación del cuadro 18.3 resulta en la siguiente ecuación cúbica para cadaintervalo:

r / \ f¡"(Xi~l) , x 3 , fi (Xi) , N 3

fi(x) = 7 - 1 Ax¡ - xY + — -(x - Xi -xY

f(*i-Ó /"(*/-! )(*/ - X¡-Ó'

+

x¡ — Xi-\ 6

/(*/) f"(Xi)(Xj -Xj.Q

Xi — 6

(x¡ - x)

( A - X , - , ) ( 18 .36)

Esta ecuación con tiene sólo do s incógnitas ( las segundas derivadas al f inal de cada intervalo). Esas incógnitas se pueden evaluar mediante la siguiente ecuación:

(x ¡ - * ¡ _ I ) / V Í - Í ) + 2 ( x , + 1 - x , _ , ) / U ) + (JCi+ i - x,)f"(x¡+l

)

6 -/(*,•)]+ _6 \f(x¡-ú - f{x¡)] ( 1 8 . 3 7 )

X¡ + \ — X¡ X¡ — Xj-\



5 3 2I NTERPOLACI ÓN

C u a d r o 1 8 . 3 D E R I V A C I Ó N D E S E G M E N T A R I A S C Ú B I C A S

]a derivación (Cheney y Kincaid, 1985) se1.1 primer paso er>^ j o n d e c o m o c a d a p a r d e n o d o s s e c o t l e c t a

biisii en la observa e g u n d a derivada dentro de cada intervalo espor una cúbica, la ¡> c u a c i ó n ( 1 8 3 5 ) s e p u e d e diferenciar dosuna linea recta. L*» e s t a o b s e r v a c i ó n . C o n e s t a b a s 6 ; ) a s e g u n d avoces para verifica p r e s e n t a r con una interpolación polinomialderivada se puede í Q rden [v éase e c u a d ó n (,g 22)].de Lagrange de prv

x — x¡f¡\x) = f!Xx¡^>Xi-l ~Xi

+ fi(M)- X¡-1

Xi - x¡.(B18.3.1)

jot de la segunda derivada en cualquier pun-úonde f¡"(x) es el V i n t e r valo. Así, esta ecuación es una línealo x dentro de i-ésV g e g u n d a d e r ivada en el primer nodo/"(*,_i)recta que conecta 1 * ^£ n e , s e g u n d o n o d o f , ( x ¡ ) _con la s eg un da d e f í ^ c i ó n ( B 1 8 3 j) s e p u e d e integrar dos veces

Después, la e¿ e s ¡ o n p a r a ^ § i n e m D a r g 0 ; e s t a e x p r e .para obtener una n s t a n t e s desconocidas de integración. EsasH i ó n contendrá dos e v a l u a r a l l l a m a r a l a f u n c i ó n d e condicio-constantes se pued* ) d e b e s e r ¡ g u a l a f ( x ¡ t ) e n ^ y f ( x ) d e b e

nos de igualdad / A 1 r e a i ¡ z a r e s t a s evaluaciones, resulta laser igual &f(xt) en % b i c a :

siguiente ecuación 6

f!Xx¡fi(X) =

(x¡ - xY +f'M)

6(x ¡— Xi-

f'XxMxi-Xi-i)'

:(-v-.v,_]r

(x¡ - x)

(X - A V - l )

(B 18.3.2)

esta relación es una expresión mucho másj Ahora, es claro qu^ ^ ^ ^ ^ cúbica para el ¡-ésimo interva-

comp heada para ufl u a c i o n ( 1 8 3 5 ) s i n e m b a r g 0 , observe quelo que, digamos, la 0

, , ,

contiene solo dos "coeficientes" desconocidos, las segundas derivadas al inicio y al final del intervalo —fXx¡_,) y/"(*, )—. Deesta forma, sipodemos determinar la segunda derivada adecuada en cada nodo, la ecuación (B 18.3.2) es un polinomial de tercer orden que se puede usar para interpolar dentro del intervalo.

Las segundas derivadas se evalúan al llamar las condiciones de que las primeras derivadas en los nodos deben ser continuas:

(B18.3.3)

La ecuación (B18.3.2) puede diferenciarse con el finde dar unaexpresión para la primera derivada. Si esto se hace tanto para el(/' — l)-ésimo , como para /-ésimo intervalo y los dos resultadosse igualan de acuerdo con la ecuación (B18.3.3), resulta la siguiente relación:

(x¡ - A - ; - I ) / " ( . Y ; - I ) + 2 ( J : / + I - x¡-i)f"(x¡)

+ (Xi + \ -Xi)f"(x¡+ \)

6[f(x l+l) - f(X¡)]

Xi + \— x¡

+ IftXi-l) - f(Xi)]

Xi X¡ — \(B18.3.4)

Si la ecuación (B18.3.4) se escribe para todos los nodos interiores, n — 1 ecuaciones simultáneas resultan con n + 1 segundasderivadas desconocidas. Sin embargo, ya que ésta es unasegmentaria cúbica natural, las segundas derivadas en los nodosextremo son cero y el problema se reduce a n — 1 ecuacionescon n — 1 incógnitas. Además, observe que el sistema de ecuaciones será tridiagonal. Así, no sólo redujimos el número d e ecuaciones, sino que las forjamos en una forma que es en extremofácil de resolver (recuerde la sección 11.1.1).

Si esta ecuación es descrita para todos los nodos interiores, resulta n — 1 ecuacionessimultáneas con n — 1 incógnitas. (Recuerde que las segundas derivadas en los nodosextremo son cero.) La aplicación de esas ecuaciones se ilustra en el siguiente ejemplo.

E J E M P l ¿ ? 18.10 Segmentar ias cúbicas

E n u n c i a d o del p r o b l e m a . Ajuste por segmentarias cúbicas los mismos datos que seusaron en los ejemp los 18.8 y 18.9 (vé ase tabla 18 .1). Utilice los resultados para estimarel valor en x = 5.



Solución, lil primer paso es emplear la ecuación (18.37) para generar un conjunto de

ecuaciones simultáneas que serán utilizadas para determinar las segundas deriviiduN en

los nodos. Por ejemplo, para el primer nodo interior se usan los siguientes datos:

= 3 f(x0) = 2.5

X\ = 4.5 1

X2 = 7 f(X2) = 2.5

Estos valores pueden sustituirse en la ecuación (18.37) para dar

(4.5 - 3)/ "(3 ) + 2(7 - 3)/"( 4.5) + (7 - 4.5 )/" (7)

6 6: ( 2 . 5 - D + — - ( 2 . 5 - 1 )

1 7 - 4 . 5 4 . 5 - 3

j

| Debido a la condición de la segmentaria na tur al, /"( 3) = 0, y la ecuación se reduce a

| 8/"(4.5) + 2.5/"(7) = 9.6

| En una forma similar, la ecuación (18.37) se aplica al segundo punto interior para ob-

i tener

j 2.5/"( 4.5) + 9/"(7) = - 9 . 6

í| Estas dos ecuaciones pueden resolverse simultáneamente para dar

j /"(4.5) = 1.67909

| f " { l ) = - 1 . 5 3 3 0 8

| Estos valores se sustituyen entonces en la ecuación (18.36) , jun to con los valores de

¡ las x y las/(x), para obtenerI 1.67909 , 2.5| Mx) = (x - 3) 3 + (4.5 - x )

J 6 ( 4 . 5 - 3 ) 4 . 5 -3 '

1 1.67909(4.5 - 3)(x - 3)

_ 4 . 5 - 3 6

o

/ i ( x ) = 0.186566 (x - 3 ) 3 + 1.666667(4.5 - x ) + 0.246894 (x - 3)

Esta ecuación es la segmentaria cúbica para el primer intervalo. Sustituciones similaresse pueden hacer para desarrollar las ecuaciones para el segundo y tercer intervalo:

/ 2 ( x ) = 0.111939(7 -xf -0 .102205 ( x - 4.5) 3 - 0.299621(7 - x )

+ 1.638783(x - 4 . 5 )

y

/•, (*) =-0.1 27757 (9 - x ) 3 + 1.761027(9 - x) +().25(x - 7)




Se puede emple ar las tres ecuaciones pa ra calcular los valores dentro de cada in tervalo .Por ejemplo , el valor en x = 5 , el cual es tá dentro del segundo in terva lo , se calcula com o

/ 2 (5 ) = 0 .111939(7 - 5 ) 3 - 0 . 1 0 2 2 0 5 ( 5 - 4 . 5 ) 3 - 0 .299 621(7 - 5)

+ 1.638783(5 - 4.5) = 1.102886

Se calculan o tros valores y los resultados se graf ican e n la figura 18 .16c.

Los resultados dé lo s ejemplos 18 .8 a 18 .10 se resum en en la f igura 18 .16. Obs ervela mejora progresiva del ajuste conforme nos movemos de l ineales a cuadráticas y asegmentar ias cúbicas . También hemos sobrepuesto un polinomio de in terpolación cúbica sobre la f igura 18 .16c. Aun que la segmentar ia cúbica consis te e n una ser ie de curvasde tercer orden, la resultante d if iere de la obtenida al usar el polin om io de tercer orden.Esto se debe al hec ho de que la segmentar ia natural no requiere de segund as der ivadas enlos nodos extremo, mientras que el polinomio cúbico no t iene tal res tr icción .

1 8 . 6 . 4 A l g o r i t m o d e c ó m p u t o p a r a s e g m e n t a r i a s cú b i ca s

El método para calcular segmentar ias cúbicas , descr i to en la sección anter ior , es idealpara la imp lemen tac ión en compu tado ra . R ecuerde que , con a lgunas man ipu lac ionesin tel igentes , el método se reduce a la resolución de n — 1 ecuac iones s imu l táneas . Unbenef icio extra de la der ivación implica, com o lo especif ica la ecuació n (18 .37) , que els is tema de ecuac iones es tr id iagonal. Com o se descr ib ió en la sección 11.1, los algor i tmos

están d isponibles para resolver tales s is temas en una manera en extremo ef iciente. Laf igura 18 .18 muestra una estructura computacional que incorpora esas caracter ís t icas .

Obse rve que la ru tina en la f igura 18 .18 regresa só lo un valor in terpolado, yu , paraun valor dado de la var iable dependie nte, xu . Ésta es só lo una forma con la cual se pued eimplementar la in terpolación segmentar ia . Por ejemplo , a usted le gustar ía determinarahora los coef icientes y después realizar muchas in terpolaciones. Además, la ru tina datanto la pr imera der ivada (dy) como la segunda (dy 2 ) en xu . Aunque no es necesar iocalcular esas cantidades, son ú tiles en much as aplicaciones de la in terpolación segmentar ia .

P R O B L E M A S

18.1 Estime el logaritmo de 5 de base 10 (log 5) medianteInterpolación lineal.a) Interpole entre log 4 = 0.60206 y log 6 = 0.7781513.b) Interpole entre log4.5 = 0.6532125 y log 5.5 = 0.7403627.

Para cada una de las interpolaciones, calcule el porcentaje deerror relativo con base en el error real.IH.2 Ajuste con'un polinomio de interpolación de Newton desegundo orden para estimar log 5 por medio de los datos delproblema 18.1. Calcule el error relativo porcentual real.I N..1 Ajuste con un polinomi o de interpolación de Newton de ter-oor orden para estimar log 5 usando los datos del problema 18.1.

18.4 Dados los datos

X 1 2 2.5 3 4 5

f|x) 1 5 7 8 2 1

a) Calcule/(3.4) mediante polinomios de interpolación deNewton de orden 1 a 3. Escoja la secuencia de p untos para suestimación con el fin de obtener la mejor exactitud posible.

b) Utilice la ecuación (18.18) para estimar el error para cadapredicción.

I K.5 Dados los datos



PROBLEMA5 0 9 9

SUtíKOUTINH bpllne (x,y.n,xu,yu,dy,d2y)LOCAL en¡ fn¡ qíi rn¡ d2xn

CALL Tridiag(x,y,n,e,f,g,r)

CALL Decomp(e,f,g,n-1)

CALL 5ubst(e,f,g,r,n-1,d2x)

CALL lnterpo\(x,y,n,d2x,xu,yu,dy,d2y)

END Spline

SU&RO UTINE Tridiag (x, y, n, e, f, g, r)

f1 = 2*(xz-x 0)

01 = F X 2 - X I )

r, = 6 / (x 2-x<¡) * (yz-y :)r, = r,+6/(x-x0) * (y 0-y,;

DO 1-2, n -2

e¡ = (x¡-xH)

fi = 2*(xM-x H)

0¡ = ( x¡+i-x¡)

n = 6/(*i+i-xi)*(yM-y¡)END DO

fn-^2.*(x„-x^2)

r „_ , = 6/(x„ - xn

_,) * (y n-y n _ ,)V , = r^+6l{xnA - x n_ 2j * (y„_ 2 -y n_ ,)

END Tridiag

5UPR0U TINE Interpol (x,y,n.d?,x.xu,yu,dy.dZy)

fíag = O

1 = 1

DO

IFxu>x¡_^ANDxu<xJHEN

el = d2x-, _ /&/(x¡ - x¡ _,)

c2 = d2x/6/(xi-x¡_ 1)

c5 = (y _ /(x, -X_,) - d2x¡ _, * (x, - x¡ _ ,)/6

c4 = (y/(x¡ - x¡ _ ,) - d2x¡ * (x ¡ - x¡ _ ,)/6

t1 = d * (x,-xuf

t 2 = c 2 * ( x u - X J _ , J 3

t 3 = c3 * (x¡ - xu)1 A = c4*(xu-x¡_i)yu=t1 + t2 + t3 + t4

t í = -3*c1 *(x¡-xuf

t2 = 3*c2*(xu-x¡_1f

13 =-c3

t4 = c4

dy = t1 + t2 + t3 + t4

t í = 6* el* (x ¡ - xu)

t2 = 6*c2*(xu-xi_1)

d2y = t í + t 2fíag = í

EL5E

i = i + 1

END IFIFi = n + WK flag = 1 EXIT

EN D DO

IF flag = O THEN

PRINT "outside range"

pause

END \F

END Interpol

F I G U R A 1 8 . 1 8

Algoritmo para interpolación segmentaria cúbica.



5 3 6 INTE RPO LACIÓ N

X 1 2 3 5 6

/|x) 4.75 4 5.25 19.75 36

Calcu le f(A) usando po l inomios de in t e rpo lac ión de Newton deo rden 1 a 4 . Esco ja sus pun tos base para ob tene r una bu ena exac t i tud. ¿Qué le indican los resul tados con respecto al orden que seusó de los pol inomios para generar los datos en la tabla?

18.6 Rep i t a lo s p rob lem as 18 .1 has t a 18 .3 med ian te po l inomiosd e L ag r an g e .

18.7 Rep i t a e l p rob lema 18 .5 u sando e l po l inomio de Lag rangede orden 1 a 3 .18.8 E m p l e e i n t e r p o l ac i ó n i n v e r s a u s an d o i n t e r p o l ac i ó n d e

po l inom ios cúb ica y b i secc ión para de te rminar e l valo r dex q u eco r re sponda af(x) = 0 .3 para lo s s igu ien te s dato s t abu lados ,

X 1 . 2 3 4 5 6 7

f(x) 1 0.5 0 .3 3 3 3 0.25 0 .2 0 . 1 6 6 7 0 . 1 4 2 9

18.9 Em plee in t e rpo lac ión inve r sa para de te rminar e l valo r dexque co r re sponda af(x) = 0.93 para los s iguien tes datos tabula

dos ,

X 0 1 2 3 4 5

f(x) 0 0.5 0 .8 0.9 0 .941 17 6 0 .9 61 5 3 8

Observe que los valores en la tabla fueron generados con la fun

c i ó n / ( x ) = x2

l( \ + x2

).a) Dete rmin e en fo rma anal í t ic a e l valo r co r r ec to .

b) Use in te rpo lac ión cuad rát i ca y l a fó rmu la cuad rát i ca parade te rminar numér icamen te e l valo r .

c) Use in terpolación cúbica y bisección para determinar el valord e m an e r a n u m é r i c a .

18.10 Desar ro l l e s egmen tar i as cuad rát i cas para lo s p r imeros 5datos en e l p rob lema 18 .4 y p red iga/ (3 .4 ) y f(2.2).

18 . 11 Desar ro l l e s egm en tar i as cúb icas para lo s dato s en e l p ro b lema 18 .5 , y a) p r e d i g af(4) y7(2.5) y b) verif ique que7^(3) y

/ 3 ( 3 ) = 5 . 2 5 .

1 8 .1 2 Dete rm ine lo s coe f i c i en te s de l a parábo la que pasa a t ravés de lo s t r e s ú l t imos pun tos de l p rob lem a 18 .4 .

18.13 Dete rm ine lo s coe f i c i en te s de la ecuac ión cúb ica que pasaa t r avés de lo s p r ime ros cuat ro pun tos de l p rob lem a 18 .5 .

1 8 .1 4 Desar ro l l e , depu re y p ruebe un p rog rama de p rueba en unlenguaje de al to n ivel o un lenguaje macro de su e lección para

imp lemen tar in t e rpo lac ión de po l inomios de Newton con baseen la f igura 18.7 .

18.15 Pruebe e l p rog rama que desar ro l ló en e l p rob lema 18 .14para lo s mismo s cálcu los de l e j emplo 18 .5 .

18.16 Use e l p rog ram a que desar ro l ló en e l p rob lem a 18 .14 parare so lve r lo s p rob lem as 18 .1 has t a 18 .3 .

18.17 Use e l p rog rama que desar ro l ló en e l p rob lem a 18 .14 pararesolver los problemas 18.4 y 18.5 . En e l problema 18.4 , u t i l icetodos lo s dato s para desar ro l l ar po l inomios de p r imer o rden a lqu in to . Para ambos p rob lem as , g raf ique e l er ro r e s t imad o con t r ael orden.

18.18 Desar ro l l e , depu re y p ruebe un subp ro g rama en l enguaje s

de al to n ivel o en un lenguaje macro de su e lección para implemen tar in t e rpo lac ión de Lag range . Con base en e l p seudoc ód igode la f igura 18.11. Pruébe lo con lo s mismos dato s de l e j emplo18.7 .

18.19 Una apl icación út i l de la in terpolación de Lagrange es l lam a d a tabla de consulta. Co mo su nom bre lo implica, és ta involucra"consu l t ar " un valo r in t e rmed io de l a t ab la . Para desar ro l l ar unalgor i tmo como tal , la tabla de x y los valores def(x) se guardanpr imero en un par de ar r eg los un id imens ion ale s . Esos valo re s s epasan entonces a una función junto con e l valor dex que us t edquiera evaluar . La función entonces e jecuta dos tareas . Pr imero,va por la tabla hasta que encu entra e l in tervalo dentro del cual es tál a incógn i t a . Después se ap l i ca una t écn ica como l a in t e rpo lac iónde Lag rang e para de te rminar e l valo r adecua do def(x). D e s a r r o

l l e t a l func ión u san do un po l inom io de Lag ran ge cúb ico para r eal izar la in terpolación. Para valores in termedios es una excelenteselección, ya que la incógni ta es tará ubicada en e l in tervalo medio de los cuatro puntos necesar ios para generar la cúbica. Seacuidadoso con los in tervalos pr imero y úl t imo donde és te no es e lcaso . También t enga su cód igo para de tec tar cuándo r equ ie re e lusuar io un valor fuera del rango de lasx. Para tales casos , la func ión debe r í a mos t r ar un e r ro r de mensaje . P ruebe su p rog ramapara f(x) = ln x median te dato s desde x = 0, 1, 2 , . . . , 10.

18.20 Desar ro l l e , depu re y p ruebe un subp rog rama en un l en guaje de al to n ivel o un lenguaje macro de su e lección para implementar in terpolación segmentar ia cúbica con base en la f igura18.18. P ruebe e l p rog rama con lo s mismos dato s de l e j emplo

18.10.18 . 21 Use el sof tware que desarrol ló en e l problema 18.20 paraajustar segmentar ias cúbicas a t ravés de los datos en los problema s 18 .4 y 18 .5 . Para amb os casos , p re d ig a/ (2 .25 ) .

Chapra17-18

Documents

Transcript of Chapra17-18