Non-linear Curve Fitting

Transcribed from the printed issue; not yet reviewed. Read it in the PDF of the issue, page 126

David Appleton is right to observe (Vector Vol. 8 No. 3) that curve fitting routines for non-linear functions should be included in ASL, and indeed more generally there should be awareness that although it is only linear fits which come as a primitive with APL, the domino function nevertheless allows fitting of a much wider range of non-linear functions at the cost of a very modest amount of programming.

I take issue with David in using a technique which requires the analyst to work out his partial derivatives in algebra - surely a denial of what a computer is for, quite apart from the unnecessary wear on the execute key!

Take as a starting point David’s example function:

y=a-(a-b)exp(-cx)

and define as FN its realisation as an APL function (C stands for coefficients):

[0]   Z←C FN X
[1]   Z←C[1]+(-/C[2 1])×*-C[3]×X

Partial derivatives of a function F are estimated by defining an operator:

[0]   Z←C(F PDERIV X)N;T
[1]   ⍝ N is index of coefficient whose partial deriv. is reqd
[2]   Z←(((C+T×N=⍳⍴C)F X)-C F X)÷T←∆X3[N]

David in his derivative function uses a blanket .01% for ∆X, however intervals can be defined for each coefficient separately as the rightmost part of the above function line implies:

      ∆X3
0.01 0.01 0.001

Now use domino to do the fit:

[0]   Z←X(F FIT Y)C
[1]   Z←C+(Y-C F X)⌹⊃[1](⊂C)F PDERIV X)¨⍳⍴C

Given the following observed readings which are approximately the same as David’s:

      X
0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19
      Y
21 27 29 31 30 37 35 32 34 37 39 37 37 40 36 37 40 41 40 38

and a first guess at the coefficients:

      C0
40 20 0.2

running the function FIT gives:

      X FN FIT Y C0
39.303 22.524 0.20628

Now develop a means of finding an iterative solution with prescribed accuracy. Start by defining a completely general operator which executes a dyadic function P feeding back the result until a condition defined by the predicate function Q is attained:

[0]   Z←L(P UNTIL Q)R
[1]   →0 IF Q Z←R
[2]   Z←L P UNTIL Q L P R

The function Q defines a stopping criterion which might be that all the coefficients are within EPS of the previous iteration. The function ALLNEAR uses the P and L of UNTIL:

[0]   Z←ALLNEAR C
[1]   Z←∧/EPS>|C-C P L

      EPS
0.00001

A small amendment must be made to FIT:

[0]   Z←A(F FIT1 Y)B
[1]   Z←R+(Y-R F L)⌹⊃[1](⊂R)(F PDERIV L)¨⍳⍴R

which now has dummy arguments A and B since it is the derived function (F FIT1 Y) which is the left operand of UNTIL and the L and R are those of UNTIL. The iterated solution is:

      X FN FIT1 Y UNTIL ALLNEAR C0
39.281 22.493 0.20907

To obtain a trace of the steps towards convergence add ⎕← at an appropriate point in the operator UNTIL:

[2]   Z←L P UNTIL Q L P ⎕←R

and rerun:

      X FN FIT1 Y UNTIL ALLNEAR C0
40 20 0.2
39.303 22.524 0.20628
39.289 22.502 0.2086
39.283 22.495 0.20899
39.282 22.494 0.20905
39.281 22.493 0.20906
39.281 22.493 0.20907

The above analysis establishes the point that only a small amount of APL programming need be added to domino to generalise from linear to non-linear fitting in a manner consistent with the overall programming philosophies of APL2.