Transcript PPT

Nonlinear Regression
Civil Engineering Majors
Authors: Autar Kaw, Luke Snyder
http://numericalmethods.eng.usf.edu
Transforming Numerical Methods Education for STEM
Undergraduates
http://numericalmethods.eng.usf.edu
1
Nonlinear Regression
http://numericalmethods.eng.usf.edu
Nonlinear Regression
Some popular nonlinear regression models:
bx
1. Exponential model: ( y  ae )
( y  ax b )
ax 

y

3. Saturation growth model: 

b

x


4. Polynomial model:
( y  a0  a1x  ...  amx m )
2. Power model:
3
http://numericalmethods.eng.usf.edu
Nonlinear Regression
Given n data points ( x1, y1), ( x 2, y 2), ... , ( xn, yn ) best fit y  f (x )
to the data, where f (x ) is a nonlinear function of x .
( xn , yn )
( x2 , y2 )
( xi , yi )
(x , y )
1
y  f (x)
yi  f ( xi )
1
Figure. Nonlinear regression model for discrete y vs. x data
4
http://numericalmethods.eng.usf.edu
Regression
Exponential Model
5
http://numericalmethods.eng.usf.edu
Exponential Model
Given ( x1, y1), ( x 2, y 2), ... , ( xn, yn ) best fit y  ae bx to the data.
( xn , yn )
( x2 , y2 )
( xi , yi )
(x , y )
1
y  aebx
yi  f ( xi )
1
Figure. Exponential model of nonlinear regression for y vs. x data
6
http://numericalmethods.eng.usf.edu
Finding constants of Exponential Model
The sum of the square of the residuals is defined as
n

Sr   yi  ae
i 1

bxi 2
Differentiate with respect to a and b
n
S r
  2 y i  ae bxi  e bxi  0
a
i 1





n
S r
bxi
bxi
  2 y i  ae
 axi e  0
b
i 1
7

http://numericalmethods.eng.usf.edu
Finding constants of Exponential Model
Rewriting the equations, we obtain
n
  yi e
bxi
i 1
n
 y i xi e
i 1
8
bxi
n
 a  e 2bxi  0
i 1
n
 a  xi e
i 1
2bxi
0
http://numericalmethods.eng.usf.edu
Finding constants of Exponential Model
Solving the first equation for a yields
n
a
 yi e
i 1
n
e
bxi
2bxi
i 1
Substituting a back into the previous equation
n
n
 y i xi e
i 1
bxi

 yi e
i 1
n
e
bxi
n
 xi e
2bxi i 1
2bxi
0
i 1
Nonlinear equation in terms of b
The constant b can be found through numerical methods such
as the bisection method or secant method.
9
http://numericalmethods.eng.usf.edu
Example 1-Exponential Model
Many patients get concerned when a test involves injection of a
radioactive material. For example for scanning a gallbladder, a
few drops of Technetium-99m isotope is used. Half of the
techritium-99m would be gone in about 6 hours. It, however,
takes about 24 hours for the radiation levels to reach what we
are exposed to in day-to-day activities. Below is given the
relative intensity of radiation as a function of time.
Table. Relative intensity of radiation as a
function of time.
t(hrs)

10
0
1
3
5
7
9
1.000
0.891
0.708
0.562
0.447
0.355
http://numericalmethods.eng.usf.edu
Example 1-Exponential Model cont.
The relative intensity is related to time by the equation
t
  Ae
Find:
a) The value of the regression constants A and 
b) The half-life of Technium-99m
c) Radiation intensity after 24 hours
11
http://numericalmethods.eng.usf.edu
Plot of data
12
http://numericalmethods.eng.usf.edu
Constants of the Model
  Aet
The value of λ is found by solving the nonlinear equation
n
n
f      i t i e
ti
i 1

ti

e
 i
i 1
n
2ti
e

n
2ti
t
e
0
 i
i 1
i 1
n
A
t

e
 i
i
i 1
n
2 t i
e

i 1
13
http://numericalmethods.eng.usf.edu
Setting up the Equation in MATLAB
n
n
f      i t i e
i 1
ti

ti

e
 i
i 1
n
e
n
2ti
t
e
0
 i
2ti i 1
i 1
t (hrs)
0
1
3
5
7
9
γ
1.000 0.891 0.708 0.562 0.447 0.355
14
http://numericalmethods.eng.usf.edu
Setting up the Equation in MATLAB
n
n
f      i t i e
i 1
ti

 ie
i 1
n
ti
2ti
e

n
 ti e
i 1
2ti
0
  0.1151
i 1
t=[0 1 3 5 7 9]
gamma=[1 0.891 0.708 0.562 0.447 0.355]
syms lamda
sum1=sum(gamma.*t.*exp(lamda*t));
sum2=sum(gamma.*exp(lamda*t));
sum3=sum(exp(2*lamda*t));
sum4=sum(t.*exp(2*lamda*t));
f=sum1-sum2/sum3*sum4;
15
http://numericalmethods.eng.usf.edu
Calculating the Other Constant
The value of A can now be calculated
6
A
 e
i 1
6
ti
i
2 ti
e

 0.9998
i 1
The exponential regression model then is
  0.9998 e
16
0.1151t
http://numericalmethods.eng.usf.edu
Plot of data and regression curve
  0.9998 e0.1151t
17
http://numericalmethods.eng.usf.edu
Relative Intensity After 24 hrs
The relative intensity of radiation after 24 hours
  0.9998  e
0.115124 
 6.3160  102
6.316  102
 100  6.317%
This result implies that only
0.9998
radioactive intensity is left after 24 hours.
18
http://numericalmethods.eng.usf.edu
Homework
1.
2.
3.
19
What is the half-life of technetium
99m isotope?
Compare the constants of this
regression model with the one
where the data is transformed.
Write a program in the language
of your choice to find the
constants of the model.
http://numericalmethods.eng.usf.edu
THE END
http://numericalmethods.eng.usf.edu
20
http://numericalmethods.eng.usf.edu
Polynomial Model
Given ( x1, y1), ( x 2, y 2), ... , ( xn, yn ) best fit y  a0  a1 x  ...  am x
m
to a given data set. (m  n  2)
( xn , yn )
( x2 , y2 )
( xi , yi )
(x , y )
1
y  a0  a1 x    am x m
yi  f ( xi )
1
Figure. Polynomial model for nonlinear regression of y vs. x data
21
http://numericalmethods.eng.usf.edu
Polynomial Model cont.
The residual at each data point is given by
Ei  y i  a 0  a1 xi  . . .  a m xim
The sum of the square of the residuals then is
n
S r   Ei2
i 1
n

  y i  a 0  a1 xi  . . .  a m xim

2
i 1
22
http://numericalmethods.eng.usf.edu
Polynomial Model cont.
To find the constants of the polynomial model, we set the derivatives
with respect to ai where i  1, m, equal to zero.
n
S r
  2. yi  a0  a1 xi  . . .  am xim (1)  0
a0 i 1


n
S r
  2. yi  a0  a1 xi  . . .  am xim ( xi )  0
a1 i 1






n
S r
  2. yi  a0  a1 xi  . . .  am xim ( xim )  0
am i 1

23

http://numericalmethods.eng.usf.edu
Polynomial Model cont.
These equations in matrix form are given by


 n

 n
   xi 
  i 1 
. . .

 n m 
  xi 
 i 1 
 n

  xi 
 i 1 
 n 2
  xi 
 i 1 
. . .
. .
. .
. .
 n m1 
  xi  . .
 i 1



 n m 
.  xi  a
 i 1    0

 n m1  a1
.  xi 
 i 1
 . .
. . .  a
 m
n


.  xi2 m  
 i 1


 n

   yi
  ni 1

xi yi


.  i 1
 . . .
 
n
 xim yi
 i 1











The above equations are then solved for a0 , a1 ,, am
24
http://numericalmethods.eng.usf.edu
Example 2-Polynomial Model
Regress the thermal expansion coefficient vs. temperature data to
a second order polynomial.
25
Coefficient of
thermal
expansion, α
(in/in/oF)
80
6.47×10−6
40
6.24×10−6
−40
5.72×10−6
−120
5.09×10−6
−200
4.30×10−6
−280
3.33×10−6
−340
2.45×10−6
6.00E-06
5.00E-06
(in/in/o F)
Temperature, T
(oF)
7.00E-06
Thermal expansion coefficient, α
Table. Data points for
temperature vs α
4.00E-06
3.00E-06
2.00E-06
-400
-300
-200
1.00E-06
-100
0
100
200
Temperature, o F
Figure. Data points for thermal expansion coefficient vs
temperature.
http://numericalmethods.eng.usf.edu
Example 2-Polynomial Model cont.
We are to fit the data to the polynomial regression model
α  a0  a1T  a 2T 2
The coefficients a0 ,a1 , a2 are found by differentiating the sum of the
square of the residuals with respect to each variable and setting the
values equal to zero to obtain

 n
  n 2 
 n

n
T
T







i
i


  i 
i

1
i

1





 a0   i 1

n
n
n
n
 





a    T  
2
3
T
T
T









i
i
i

 1    i i 
i 1
i 1
i 1
i 1






 n



a
n


n
n
 2
2
 T 2   T 3   T 4  
 Ti  i 


i
i
i
 

 i 1

 i 1   i 1   i 1  
26
http://numericalmethods.eng.usf.edu
Example 2-Polynomial Model cont.
The necessary summations are as follows
Table. Data points for temperature vs.
Temperature, T
(oF)
Coefficient of
thermal expansion,
α (in/in/oF)
α
7
T
i 1
T
6.47×10−6
40
6.24×10−6
−40
5.72×10−6
−120
5.09×10−6
i 1
−200
4.30×10−6
7
3.33×10−6
−340
2.45×10−6
2.5580 105
3
  7.0472  10 7
4
 2.1363 1010
7
80
−280
2
i
i 1
i
7
T
i

i 1
i
7
T 
i 1
i
7
T
i 1
27
 3.3600  10 5
i
2
i
  2.6978  10 3
 i 8.5013  10 1
http://numericalmethods.eng.usf.edu
Example 2-Polynomial Model cont.
Using these summations, we can now calculate a0 ,a1 , a2
 7.0000

2
 8.600  10
 2.5800  10 5

 8.6000  10 2
2.5800  10 5
 7.0472  10 7
2.5800  10 5  a 0   3.3600  10 5 



 7.0472  10 7   a1    2.6978  10 3 
2.1363  1010  a 2   8.5013  10 1 
Solving the above system of simultaneous linear equations we have
6
a 0   6.0217  10 
 a    6.2782  10 9 

 1 

11
a 2   1.2218  10 
The polynomial regression model is then
α  a0  a1T  a 2T 2
 6.0217  10 6  6.2782  10 9 T  1.2218  10 11 T 2
28
http://numericalmethods.eng.usf.edu
Linearization of Data
To find the constants of many nonlinear models, it results in solving
simultaneous nonlinear equations. For mathematical convenience,
some of the data for such models can be linearized. For example, the
data for an exponential model can be linearized.
As shown in the previous example, many chemical and physical processes
are governed by the equation,
y  aebx
Taking the natural log of both sides yields,
ln y  ln a  bx
Let z  ln y and a 0  ln a
We now have a linear regression model where z  a0  a1 x
(implying)
29
a  eao with a1  b
http://numericalmethods.eng.usf.edu
Linearization of data cont.
Using linear model regression methods,
a1 
n
n
n
i 1
i 1
i 1
n  xi z i   xi  z i


n xi2    xi 
i 1
 i 1 
n
_
n
2
_
a 0  z  a1 x
Once ao , a1 are found, the original constants of the model are found as
b  a1
a  e a0
30
http://numericalmethods.eng.usf.edu
Example 3-Linearization of data
Many patients get concerned when a test involves injection of a radioactive
material. For example for scanning a gallbladder, a few drops of Technetium99m isotope is used. Half of the technetium-99m would be gone in about 6
hours. It, however, takes about 24 hours for the radiation levels to reach what
we are exposed to in day-to-day activities. Below is given the relative intensity
of radiation as a function of time.
t(hrs)

0
1
3
5
7
9
1.000
0.891
0.708
0.562
0.447
0.355
1
Relative intensity of radiation, γ
Table. Relative intensity of radiation as a function
of time
0.5
0
0
5
Time t, (hours)
10
Figure. Data points of relative radiation intensity
vs. time
31
http://numericalmethods.eng.usf.edu
Example 3-Linearization of data cont.
Find:
a) The value of the regression constants A and 
b) The half-life of Technium-99m
c) Radiation intensity after 24 hours
The relative intensity is related to time by the equation
  Aet
32
http://numericalmethods.eng.usf.edu
Example 3-Linearization of data cont.
Exponential model given as,
  Ae t
ln    ln  A  t
Assuming z  ln  , ao  ln  A and a1   we obtain
z  a0  a1t
This is a linear relationship between z and t
33
http://numericalmethods.eng.usf.edu
Example 3-Linearization of data cont.
Using this linear relationship, we can calculate a0 , a1
a1 
n
n
n
i 1
i 1
i 1
where
n  t i zi   t i  zi
 n 
2
n t1    ti 
i 1
 i 1 
n
and
2
a0  z  a1t
a
1
a0
Ae
34
http://numericalmethods.eng.usf.edu
Example 3-Linearization of Data cont.
Summations for data linearization are as follows
Table. Summation data for linearization of data model
ti
i
1
2
3
4
5
6
0
1
3
5
7
9

25.000
i
zi  ln  i
1
0.891
0.708
0.562
0.447
0.355
0.00000
−0.11541
−0.34531
−0.57625
−0.80520
−1.0356
0.0000
−0.11541
−1.0359
−2.8813
−5.6364
−9.3207
0.0000
1.0000
9.0000
25.000
49.000
81.000
−2.8778
−18.990
165.00
ti zi
t
2
i
With n  6
6
t
i 1
6
z
i 1
 2.8778
i
6
t z
i 1
i i
6
t
i 1
35
 25.000
i
2
i
 18.990
 165.00
http://numericalmethods.eng.usf.edu
Example 3-Linearization of Data cont.
Calculating a0 , a1
a1 
6 18.990  25 2.8778
2
6165.00  25
 0.11505
a0 
 2.8778
25
  0.11505
6
6
 2.6150  10 4
Since
a0  ln  A
A  e a0
e
2.6150104
 0.99974
also
  a1  0.11505
36
http://numericalmethods.eng.usf.edu
Example 3-Linearization of Data cont.
0.11505t
Resulting model is   0.99974  e
1
  0.99974  e 0.11505t
Relative
Intensity
0.5
of
Radiation,
0
0
5
10
Time, t (hrs)
Figure. Relative intensity of radiation as a function of
temperature using linearization of data model.
37
http://numericalmethods.eng.usf.edu
Example 3-Linearization of Data cont.
The regression formula is then
  0.99974  e 0.11505t
1



b) Half life of Technetium 99 is when
2
1
0.99974  e 0 .11505t  0.99974e 0 .115050 
2
e 0 .11508t  0.5
 0.11505t  ln 0.5
t  6.0248 hours
38
t 0
http://numericalmethods.eng.usf.edu
Example 3-Linearization of Data cont.
c) The relative intensity of radiation after 24 hours is then
  0.99974e 0.1150524
 0.063200
6.3200 10 2
100  6.3216% of the radioactive
This implies that only
0.99983
material is left after 24 hours.
39
http://numericalmethods.eng.usf.edu
Comparison
Comparison of exponential model with and without data linearization:
Table. Comparison for exponential model with and without data
linearization.
With data linearization
(Example 3)
Without data linearization
(Example 1)
A
0.99974
0.99983
λ
−0.11505
−0.11508
Half-Life (hrs)
6.0248
6.0232
Relative intensity
after 24 hrs.
6.3200×10−2
6.3160×10−2
The values are very similar so data linearization was suitable to
find the constants of the nonlinear exponential model in this
case.
40
http://numericalmethods.eng.usf.edu
Additional Resources
For all resources on this topic such as digital audiovisual
lectures, primers, textbook chapters, multiple-choice
tests, worksheets in MATLAB, MATHEMATICA, MathCad
and MAPLE, blogs, related physical problems, please
visit
http://numericalmethods.eng.usf.edu/topics/nonlinear_r
egression.html
THE END
http://numericalmethods.eng.usf.edu