Curve fitting--4 S. Arlinghaus Logistic curve fitting. The exponential function leads to a situation of unbounded population growth. Assumption for the exponential: The rate of population growth or decay at any time t is proportional to the size of the population at t. Symbolically: let Y_t [N from last time] represent the size of a population at time t. The rate of growth of Y_t is proportional to Y_t: dY_t / dt = kY_t where k is a constant of proportionality. This assumption led to the curve Y_t = Y_(t0) e^(kt) ---- a situation in which growth is unbounded as t becomes large. Bounded growth--the logistic function. Assume now that in reality, when the population gets large, environmental factors dampen growth --- a population-environment dynamic. To express this idea symbolically: assume that eventually the growth rate decreases ---- dY_t / dt eventually decreases. assume that population size is limited to some maximum, q, where 0 < Y_t < q. As Y_t approaches q, dY_t / dt approaches 0, so that population size tends to become stable as t approaches infinity. The result is a model that is exponential in shape initially and includes effects of environmental resistance in larger populations. A typical geometric characterization is shown below. t The various coordinates in this figure are derived from expressing the assumptions of bounded growth algebraically. One such expression is: dY_t / dt = kY_t * (q-Y_t)/q This equation meets the conditions of bounded growth because in the factor (q-Y_t)/q when Y_t is small (q-Y_t)/q is close to 1 and the growth is therefore close to the exponential in the early stages; when Y_t is large--that is, close to q-- (q-Y_t)/q is close to 0, and the growth rate dy_t / dt tapers off to the horizontal asymptote. The factor (q-Y_t)/q, which is what alters the symbolic form of the logistic from the exponential, acts as a damper to growth--as required. ************************************************************************************************************* To understand where the various quatities in the graph above come from, consider the following set of manipulations. Replace k/q by K so that the logistic equation now reads: dY_t / dt = KY_t(q-Y_t) in which the rate of growth is proportional to the product of the population size and the difference between the maximum size and the population size. Solve this latter equation for Y_t ---- separate the variables of the differential equation (dY_t)/(Y_t(q-Y_t)) = K dt to yield INT (dY_t)/(Y_t(q-Y_t)) = INT K dt. Use a Table of INtegrals on the rational form in the left-hand integral: 1/q * ln|(Y_t)/(q-Y_t)| = Kt +C or, ln|Y_t/(q-Y_t)| = qKt +qC. Because Y_t > 0 and q-Y_t > 0, ln (Y_t/(q-Y_t)) = qKt +qC. Therefore, Y_t/(q-Y_t) = e^(qKt)e^(qC). Replace e^qC by A. Therefore Y_t/(q-Y_t) = Ae^(qKt) Y_t=(q-Y_t)Ae^(qKt) Y_t=qAe^(qKt)-Y_tAe^(qKt) Y_t+Y_tAe(qKt)=qAe^(qKt) Y_t = (qAe^(qKt)/(Ae^(qKt) + 1)) Now divide top and bottom of the last equation by Ae^(qKt) [equivalent to multiplying by 1], so that Y_t = q/(1+(1/(Ae^(qKt)))) = q/(1+ 1/A * e^(-qKt)). Replace 1/A by a and -qK by b producing a common form for the logistic function: Y_t = q/ (1 + ae^(bt)) with b < 0 because b = -qK, and q, K > 0. ********************************************************************************************************* Facts about the graph of the logistic equation. a. The line Y_t = q is a horizontal asymptote for the graph. This is so because, for b<0, lim (t --> inf) q/(1+ae^(bt)) ------> q/(1+a(0)) = q. Can the logistic curve cross this horizontal asymptote to produce one or multiple crossings? That is, can the curve approach the asymptote from above, or might the curve oscillate back and forth across the asymptote while closing down on it? This idea can be expressed in notation by setting Y_t = q in the logistic equation and solving to see if there are any points that lie simultaneously on the curve and the asymptote. Can it be that: Y_t = Y_t / (1+ae^(bt))? Or, 1=1+ae^(bt)? Or, ae^(bt) =0--or that a=0--no, becasue a=1/A. Or that e^(bt)=0--no. Thus, the logistic growth curve cannot cross the horizontal asymptote so that it approaches it entirely from one side--in this case, from below. b. Find the coordinates of the inflection point of the logistic curve. The inflection point is the point at which the curve changes from concave up to concave down. It is also the point at which the growth rate is a maximum--to find this maximum differentiate one form of the logistic equation: logistic equation: dY_t /dt = KY_t(q-Y_t)=KqY_t - K(Y_t)^2 derivative: d^2Y_t /dt^2 = Kq - 2KY_t Set this last equation equal to zero. Kq - 2KY_t = 0. Therefore, Y_t = q/2. This is the vertical coordinate of the inflection point of the curve for Y_t, the logistic curve -- dY_t /dt is increasing to the left of q/2 and increasing to the right (from evidence of second derivatives). The value found is a maximum. Horizontal component: To find t, put Y_t = q/2 in the logistic equation and solve: q/2 = q/(1 + ae^(bt)) Solving, 1 + ae^(bt) = 2; e^(bt)=1/a; e^(-bt)=a; -bt=ln a, and t = (ln a) / (-b) Thus, the labelling of the logistic curve is as given earlier-- the height of the inflection point is half the height of the horizontal asymptote. LOGISTIC CURVE FIT USING ACTUAL DATA. EQUATION: y=q/(1+ae^(bt)), b<0. FIND: from the given data, values for a and b. Then, use the equation on a number of t-values to graph the equation. TO FIND a and b: All that is required is two values of t, an upper bound q, and the assumption that the distribution is to be logistic. EXAMPLE: Inland coal consumption in India--'000 tons Data of K. Hornbarger (545 student in 1992). t y least sqs year 0 is exponential fit 1970; year 16 is 1987 projected beyond 1987 0 77435 74332.97938 1 81483 78484.54458 2 82077 82867.97851 3 87569 87496.23125 4 92759 92382.97616 5 102426 97542.65027 6 105436 102990.4969 7 101409 108742.6108 8 106165 114815.9855 9 108040 121228.5638 10 127545 127999.2905 11 136160 135148.1685 12 145569 142696.3179 13 148470 150666.0384 14 163950 159080.8751 15 174324 167965.6881 16 184790 177346.7261 Now suppose instead that estimates based on remaining reserves require that inland consumption never exceed 250000 tons of coal per year. Suppose also that the last year for which we have actual data is the year 1987. Assume a logistic representation in the form y=250000/(1+ae^(bt)) and find a and b. In 1970, when t=0, y=77435 Thus, 77435=250000/(1+a) So that 1+a=250,000/77435= Thus, a=2.228514 Now, find b. Use other actual data endpoint from 1987. In 1987, when t=16, y=184790 184790=250000/(1+2.228514*exp(16*b)) 1+2.228514*exp(16*b)=250,000/184,790= 1.352887061 exp(16*b)=(1.352887-1)/2.228514= 0.158350811 16b=ln(0.158350)= 0.459637587 b=(ln(0.158350))/16= -0.115184219 Logistic equation fit using actual data: y=250,000/(1+2.228514e^(-0.11518*t)) t logistic exponential actual 0 77435.0057 74332.97938 77435 1 83722.23968 78484.54458 81483 2 90252.97004 82867.97851 82077 3 96995.96015 87496.23125 87569 4 103915.0238 92382.97616 92759 5 110969.682 97542.65027 102426 6 118116.0349 102990.4969 105436 7 125307.8122 108742.6108 101409 8 132497.5521 114815.9855 106165 9 139637.8469 121228.5638 108040 10 146682.5842 127999.2905 127545 11 153588.1177 135148.1685 136160 12 160314.3044 142696.3179 145569 13 166825.359 150666.0384 148470 14 173090.491 159080.8751 163950 15 179084.3105 167965.6881 174324 16 184787.0014 177346.7261 184790 17 190184.2805 187251.7037 18 195267.1697 197709.8834 19 200031.6182 208752.1621 20 204478.0126 220411.1622 21 208610.6158 232721.3281 22 212436.9698 245719.0281 23 215967.2948 259442.6616 24 219213.9089 273932.7727 25 222190.6871 289232.1697 26 224912.574 305386.0522 27 227395.1549 322442.1438 28 229654.2899 340450.8338 29 231705.8104 359465.3257 30 233565.2733 379541.7944 31 235247.7695 400739.5524 32 236767.7792 423121.2246 33 238139.0698 446752.9338 34 239374.6272 471704.4957 35 240486.618 498049.6252 36 241486.3731 525866.1544 37 242384.3918 555236.2623 38 243190.3583 586246.7177 39 243913.1707 618989.1356 40 244560.9758 653560.2476 41 245141.2101 690062.188 42 245660.6434 728602.7952 43 246125.4239 769295.9308 44 246541.1236 812261.8155 45 246912.7826 857627.3842 46 247244.9529 905526.6617 47 247541.739 956101.1578 48 247806.8375 1009500.286 49 248043.5738 1065881.805 50 248254.9361 1125412.283 51 248443.6071 1188267.593 52 248611.9937 1254633.43 53 248762.2532 1324705.86 54 248896.3185 1398691.899 Suppose instead, that the 1987 data was suspect and that you felt more confident using the 1980 data. Still, a=2.228514 But, in 1980--when t=9, y=108040. So, 108,040=250,000/(1+2.228514*e^(9b)) (1+2.228512*e^(9b))=250,000/108,040= 2.313957793 Thus, e^(9b)=(2.313957-1)/2.228512= 0.589611813 So, 9b=ln(0.589611) -0.528292282 b= -0.058699142 So, the logistic equation in this case is: y=250,000/(1+2.228514e^(-0.05869914*t)) t logistic 2 logistic 1 exponential actual 0 77435.0057 77435.0057 74332.97938 77435 1 80607.00596 83722.23968 78484.54458 81483 2 83845.8089 90252.97004 82867.97851 82077 3 87147.79645 96995.96015 87496.23125 87569 4 90508.98645 103915.0238 92382.97616 92759 5 93925.04538 110969.682 97542.65027 102426 6 97391.3053 118116.0349 102990.4969 105436 7 100902.7849 125307.8122 108742.6108 101409 8 104454.2146 132497.5521 114815.9855 106165 9 108040.0652 139637.8469 121228.5638 108040 10 111654.5804 146682.5842 127999.2905 127545 11 115291.8122 153588.1177 135148.1685 136160 12 118945.6588 160314.3044 142696.3179 145569 13 122609.9054 166825.359 150666.0384 148470 14 126278.2657 173090.491 159080.8751 163950 15 129944.4254 179084.3105 167965.6881 174324 16 133602.085 184787.0014 177346.7261 184790 17 137245.0036 190184.2805 187251.7037 18 140867.0408 195267.1697 197709.8834 19 144462.1979 200031.6182 208752.1621 20 148024.6565 204478.0126 220411.1622 21 151548.8151 208610.6158 232721.3281 22 155029.3221 212436.9698 245719.0281 23 158461.1061 215967.2948 259442.6616 24 161839.4017 219213.9089 273932.7727 25 165159.7721 222190.6871 289232.1697 26 168418.1269 224912.574 305386.0522 27 171610.7368 227395.1549 322442.1438 28 174734.2426 229654.2899 340450.8338 29 177785.662 231705.8104 359465.3257 30 180762.3905 233565.2733 379541.7944 31 183662.2007 235247.7695 400739.5524 32 186483.2364 236767.7792 423121.2246 33 189224.0046 238139.0698 446752.9338 34 191883.3646 239374.6272 471704.4957 35 194460.5141 240486.618 498049.6252 36 196954.974 241486.3731 525866.1544 37 199366.571 242384.3918 555236.2623 38 201695.4187 243190.3583 586246.7177 39 203941.8984 243913.1707 618989.1356 40 206106.6384 244560.9758 653560.2476 41 208190.4932 245141.2101 690062.188 42 210194.5229 245660.6434 728602.7952 43 212119.9724 246125.4239 769295.9308 44 213968.2506 246541.1236 812261.8155 45 215740.9111 246912.7826 857627.3842 46 217439.6328 247244.9529 905526.6617 47 219066.2015 247541.739 956101.1578 48 220622.4923 247806.8375 1009500.286 49 222110.4537 248043.5738 1065881.805 50 223532.0918 248254.9361 1125412.283 51 224889.4562 248443.6071 1188267.593 52 226184.6266 248611.9937 1254633.43 53 227419.7012 248762.2532 1324705.86 54 228596.7854 248896.3185 1398691.899