RosettaCodeData/Task/Multiple-regression/Phix/multiple-regression.phix

55 lines
1.4 KiB
Plaintext

constant N = 15, M=3
sequence x = {1.47,1.50,1.52,1.55,1.57,
1.60,1.63,1.65,1.68,1.70,
1.73,1.75,1.78,1.80,1.83},
y = {52.21,53.12,54.48,55.84,57.20,
58.57,59.93,61.29,63.11,64.47,
66.28,68.10,69.92,72.19,74.46},
s = repeat(0,N),
t = repeat(0,N),
a = repeat(repeat(0,M+1),M)
for k=1 to 2*M do
for i=1 to N do
s[k] += power(x[i],k-1)
if k<=M then t[k] += y[i]*power(x[i],k-1) end if
end for
end for
-- build linear system
for row=1 to M do
for col=1 to M do
a[row,col] = s[row+col-1]
end for
a[row,M+1] = t[row]
end for
puts(1,"Linear system coefficents:\n")
pp(a,{pp_Nest,1,pp_IntFmt,"%7.1f",pp_FltFmt,"%7.1f"})
for j=1 to M do
integer i = j
while a[i,j]=0 do i += 1 end while
if i=M+1 then
?"SINGULAR MATRIX !"
?9/0
end if
for k=1 to M+1 do
{a[j,k],a[i,k]} = {a[i,k],a[j,k]}
end for
atom Y = 1/a[j,j]
a[j] = sq_mul(a[j],Y)
for i=1 to M do
if i<>j then
Y=-a[i,j]
for k=1 to M+1 do
a[i,k] += Y*a[j,k]
end for
end if
end for
end for
puts(1,"Solutions:\n")
?columnize(a,M+1)[1]