#3282 lsquares and lists of list of data to be fitted on

[Gunter Königsmann]

I have a list of data sets that are huge => running lsquares_estimates() on them to fit them on an equation would need too big amounts of time.

Normally lsquares_estimates_approximate() helps in this case as it skips the step that tries to find an exact solution. But as I have not a single data set, but a list of data sets that somehow backfires:

(%i34)  load("lsquares");
    data:[
        matrix([0,0],[1,1],[2,2]),
        matrix([2,2],[1,1],[2,2])
    ];
    eqtn:y=a*x^2+b*x+c;
    lsquares_mse(data[1],[x,y],eqtn);
    lsquares_estimates_approximate(%,[a,b,c]);
(%o30)  "/usr/local/share/maxima/branch_5_39_base_66_gbbb452f/share/lsquares/lsquares.mac"
(data)  [matrix(
        [0, 0],
        [1, 1],
        [2, 2]
    ),matrix(
        [2, 2],
        [1, 1],
        [2, 2]
    )]
(eqtn)  y=a*x^2+b*x+c
(%o33)  sum((m6[i,2]-a*m6[i,1]^2-b*m6[i,1]-c)^2,i,1,1)
    Maxima encountered a Lisp error:
     The value
       ((MTIMES SIMP) -2.0
        ((MEXPT SIMP) ((MQAPPLY SIMP ARRAY) (($DATA SIMP ARRAY) 1) 1 1) 2)
        ((MPLUS SIMP) -1.0
         ((MTIMES SIMP) -1.0 ((MQAPPLY SIMP ARRAY) (($DATA SIMP ARRAY) 1) 1 1))
         ((MTIMES SIMP) -1.0
          ((MEXPT SIMP) ((MQAPPLY SIMP ARRAY) (($DATA SIMP ARRAY) 1) 1 1) 2))
         ((MQAPPLY SIMP ARRAY) (($DATA SIMP ARRAY) 1) 1 2)))
     is not of type
       DOUBLE-FLOAT.
    Automatically continuing.
To enable the Lisp debugger set *debugger-hook* to nil.

Only running lsquares_estimates on a single entry of the list works:

(%i29)  load("lsquares");
    data:matrix([0,0],[1,1],[2,2]);
    eqtn:y=a*x^2+b*x+c;
    lsquares_mse(data,[x,y],eqtn);
    lsquares_estimates_approximate(%,[a,b,c]);
(%o25)  "/usr/local/share/maxima/branch_5_39_base_66_gbbb452f/share/lsquares/lsquares.mac"
(data)  matrix(
        [0, 0],
        [1, 1],
        [2, 2]
    )
(eqtn)  y=a*x^2+b*x+c
(%o28)  sum((data[i,2]-a*data[i,1]^2-b*data[i,1]-c)^2,i,1,3)/3
    *************************************************
      N=    3   NUMBER OF CORRECTIONS=25
           INITIAL VALUES
     F=  1.000000000000000D+01   GNORM=  1.753726191728787D+01
    *************************************************
       I  NFN     FUNC                    GNORM                   STEPLENGTH
       1    2     2.796937667198332D-01   1.961076084276657D+00   5.702144409522791D-02  
       2    3     1.464269617620298D-01   4.983306664572200D-01   1.000000000000000D+00  

---snip---

      13   14     3.397773998012277D-06   4.742450181229817D-03   1.000000000000000D+00  
      14   15     2.643557109421608D-08   2.597101638571121D-04   1.000000000000000D+00  
     THE MINIMIZATION TERMINATED WITHOUT DETECTING ERRORS.
     IFLAG = 0
(%o29)  [[a=-3.405898581949557*10^-4,b=1.000676542954337,c=-1.342533507359797*10^-4]]

I have been able to ship around this problem by putting the data set I am working on into a temporary variable before running lsquares:

load("lsquares");
data:[
    matrix([0,0],[1,1],[2,2]),
    matrix([2,2],[1,1],[2,2])
];
eqtn:y=a*x^2+b*x+c;
makelist(block([tmp],
        tmp:data[i],
        lsquares_estimates_approximate(lsquares_mse(tmp,[x,y],eqtn),[a,b,c])
    ),
    i,1,length(data)
);

...but I am convinced that the fact that I needed this is a bug.

[Robert Dodier]

Yes, this is a bug in lsquares_mse. The problem is that it introduces a temporary variable and tries to assign the data to the temporary variable, but it is incorrect. I think I see how to fix it.

[Robert Dodier]

Fixed by commit b58a547. Thanks for the pointing it out. Closing this report as fixed.

[peter moueza]

Hi, still not fixed.
I run example from Fitting equations to real measurement data in Help of wxMaxima 20.12.1 Linux Ubuntu.
When :
lsquares_estimates_approximate(
lsquares_mse(
data1,[x,y],approach1
),
[a,b,c],
initial=[0,0,0]
);
params1:%[1];

It says :
Maxima encountered a Lisp error:
Condition in / [or a callee]: INTERNAL-SIMPLE-ERROR: Zero divisor.
Automatically continuing.
To enable the Lisp debugger set debugger-hook to nil.

So in Maxima :
read_matrix;fun1(x):=3x^2-2x+7;time1:makelist(i,i,-10,10,.02)$time1:makelist(i,i,-10,10,.02)$
data1:transpose(
matrix(
time1,
makelist(fun1(i)+random(32.0)-16,i,time1)
)
)$load("lsquares")$
lsquares_estimates_approximate(
lsquares_mse(
data1,[x,y],approach1
),
[a,b,c],
initial=[0,0,0]
);

[peter moueza]

Sorry, it is OK ( I might not have evaluated all steps)