Subversion Repositories wimsdev

!goto $wims_read_parm

:def
title=2D linear regression
synonyme=linear regression, least square, line fitting
input=data2d
!exit

:proc
cnt=!linecnt $fml
pari_header=data=[$fml2];a=sum(t=1,$cnt,data[t,1]^2);\
        b=sum(t=1,$cnt,data[t,1]);\
        c=sum(t=1,$cnt,data[t,1]*data[t,2]);\
        d=sum(t=1,$cnt,data[t,2]);\
        sol=matsolve([a,b;b,$cnt],[c;d]);
result=!exec pari print(sol[1,1])\
        print(sol[2,1])

test=!linecnt $result
!if $test<2
 result=
 !exit
!endif
!distribute lines $result into a,b
result=!rawmath $a*x + $b

l1=!line 1 of $fml
!distribute item $l1 into x1,y1
xmin=$[$x1]
xmax=$[$x1]
ymin=$[$y1]
ymax=$[$y1]
insplot_data=$[$x1] $[$y1]
!for i=2 to $cnt
 l=!line $i of $fml
 !distribute item $l into x_,y_
 xmin=$[min($xmin,$x_)]
 xmax=$[max($xmax,$x_)]
 ymin=$[min($ymin,$y_)]
 ymax=$[max($ymax,$y_)]
 insplot_data=!append line $[$x_] $[$y_] to $insplot_data
!next i

xrange=$[max((abs($xmin)+abs($xmax))/10000,0.000000001)]
yrange=$[max((abs($ymin)+abs($ymax))/10000,0.000000001)]
xdiff=$[max($xrange,0.1*($xmax-($xmin)))]
x1=$[$xmin-$xdiff]
x2=$[$xmax+$xdiff]
ydiff=$[max($yrange,0.1*($ymax-($ymin)))]
y1=$[$ymin-$ydiff]
y2=$[$ymax+$ydiff]

insplot_set=size $[$picx/500], $[$picy/400]
!insplot [$x1:$x2] [$y1:$y2] 'insplot_data' notitle with points pt 4, $result notitle
ins=$ins_out
!if getins notin $ins
 result=
!endif

!exit

:output
Line equation:
!htmlmath y = $result
<p><center>
$ins
</center> <p>
Plotted data: $cnt pairs { $fml3 }
<p>
!exit