Cubic regression
- zxretrosoft
- Enthusiast
- Posts: 169
- Joined: Wed May 15, 2013 8:26 am
- Location: Czech Republic, Prague
- Contact:
Cubic regression
Hello friends,
I have a question. I need to program a higher regression than quadratic, namely, cubic regression.
To do this, we need to know how the coefficients A, B, C, D are calculated.
I found coefficients A, B for linear:
https://www.easycalculation.com/statist ... ession.php
I found coefficients A, B, C for quadratic:
https://www.easycalculation.com/statist ... ession.php
But the coefficients A, B, C, D for cubic regression I can't find anywhere
Please, do you have someone tip? How are these coefficients calculated?
Thanks in advance!
I have a question. I need to program a higher regression than quadratic, namely, cubic regression.
To do this, we need to know how the coefficients A, B, C, D are calculated.
I found coefficients A, B for linear:
https://www.easycalculation.com/statist ... ession.php
I found coefficients A, B, C for quadratic:
https://www.easycalculation.com/statist ... ession.php
But the coefficients A, B, C, D for cubic regression I can't find anywhere
Please, do you have someone tip? How are these coefficients calculated?
Thanks in advance!
Re: Cubic regression
The general case is written here:
Polynomial regression - Matrix form and calculation of estimates
I think you can find matrix calculation procedures here in the forum.
Polynomial regression - Matrix form and calculation of estimates
I think you can find matrix calculation procedures here in the forum.
PB 6.01 ― Win 10, 21H2 ― Ryzen 9 3900X, 32 GB ― NVIDIA GeForce RTX 3080 ― Vivaldi 6.0 ― www.unionbytes.de
Lizard - Script language for symbolic calculations and more ― Typeface - Sprite-based font include/module
Lizard - Script language for symbolic calculations and more ― Typeface - Sprite-based font include/module
- zxretrosoft
- Enthusiast
- Posts: 169
- Joined: Wed May 15, 2013 8:26 am
- Location: Czech Republic, Prague
- Contact:
Re: Cubic regression
Yes, thank you.
But,
it's very complicated for me. I can't deduce it. I'm just looking for 4 equations for
A = ...
B = ...
C = ...
D = ...
Nothing more, nothing less.
Exactly what is here (A=... B=... C=...) for quadratic regression:
https://www.easycalculation.com/statist ... ession.php
But,
it's very complicated for me. I can't deduce it. I'm just looking for 4 equations for
A = ...
B = ...
C = ...
D = ...
Nothing more, nothing less.
Exactly what is here (A=... B=... C=...) for quadratic regression:
https://www.easycalculation.com/statist ... ession.php
Re: Cubic regression
Here is a updated version from "alter Mann":
https://www.purebasic.fr/german/viewtop ... roximation
You can use Num_ApproxPoly (n.l,m.l,Array x.d(1),Array y.d(1), Array w.d(1),Array c.d(1)) where:
n is the number of points
m is the degree of the polynom
x is the array of x-values
y is the array of y-values
w is an array for weighting (typical 1.0)
c is the resulting array of you coefficients (A,B,C,...)
https://www.purebasic.fr/german/viewtop ... roximation
You can use Num_ApproxPoly (n.l,m.l,Array x.d(1),Array y.d(1), Array w.d(1),Array c.d(1)) where:
n is the number of points
m is the degree of the polynom
x is the array of x-values
y is the array of y-values
w is an array for weighting (typical 1.0)
c is the resulting array of you coefficients (A,B,C,...)
Code: Select all
; Author: alter Mann , adaptiert von Engeln-Müllges/Uhlig "Numerical Algorithms with C"
; Date: 23. 10. 2008
; OS: Windows
; Demo: Ja
;
; Polynomapproximation (Ausgleichpolynom)
; ( wird leicht instabil bei hohen Polynomgraden und großen Stützpunkt-Werten)
#DBL_EPSILON =2.2204460492503131e-016 ; kleinster Wert, so dass 1.0+DBL_EPSILON != 1.0
#EPS_QUAD =(#DBL_EPSILON*#DBL_EPSILON)
EnableExplicit
Procedure.d MathSqrtD (Wert.d)
Protected Wert1.d
If Wert <= 0.0
Wert1 = 0.0
Else
Wert1 = Sqr(Wert)
EndIf
ProcedureReturn Wert1
EndProcedure
Macro SQRT(wert)
MathSqrtD(wert)
EndMacro
Macro SQR2(wert)
(wert)*(wert)
EndMacro
;************************************************************************
;* ein Polynom P in der Darstellung *
;* P(x) = adA[0] + adA[1] * dX + adA[2] * dX^2 + ... + adA[n] * dX^n *
;* nach dem Hornerschema auswerten *
;* *
;* Eingabeparameter: *
;* ================= *
;* lN: Grad des Polynoms *
;* adA: [0..n]-Vektor mit den Koeffizienten des Polynoms *
;* dX: Stelle, an der das Polynom auszuwerten ist *
;* *
;* Funktionswert: *
;* ============== *
;* P(x) *
;************************************************************************
Procedure.d Num_Horner (lN.l, Array adA.d(1), dX.d) ; Hornerschema zur Polynomauswertung .............*/
Protected ldErg.d = adA(lN)
Protected i.l
For i=lN-1 To 0 Step -1
ldErg = ldErg * dX + adA(i)
Next i
ProcedureReturn ldErg
EndProcedure
;*====================================================================*
;* *
;* Num_ChoCec berechnet die Zerlegungsmatrix einer symmetrischen, *
;* positiv definiten Matrix. Die Zerlegungsmatrix wird auf mat *
;* ueberspeichert. *
;* *
;*====================================================================*
;* *
;* Eingabeparameter: *
;* ================ *
;* n Long n ( n > 0 ) *
;* Dimension von mat, *
;* Anzahl der Komponenten des b-Vektors. *
;* mat Double mat(n,n) *
;* Matrix des Gleichungssystems. *
;* *
;* Ausgabeparameter: *
;* ================ *
;* mat Double mat(n,n) *
;* Dekompositionsmatrix, die die Zerlegung von *
;* mat enthaelt. *
;* *
;* Rueckgabewert: *
;* ============= *
;* = 0 alles ok *
;* = 1 n < 1 gewaehlt *
;* = 2 Matrix ist nicht positiv definit *
;* *
;*====================================================================*
Procedure.l Num_ChoDec (n.l, Array mat.d(2))
Protected j.l, k.l, i.l
Protected sum.d
If n < 1
ProcedureReturn 1 ; n < 1 Fehler !
EndIf
If mat(0,0) < #EPS_QUAD
ProcedureReturn 2 ; mat ist nicht positiv definit !
EndIf
mat(0,0) = SQRT(mat(0,0))
For j = 1 To n-1 Step 1
mat(j,0) / mat(0,0)
Next j
For i = 1 To n-1 Step 1
sum = mat(i,i)
For j = 0 To i-1 Step 1
sum - SQR2(mat(i,j))
Next j
If sum < #EPS_QUAD
ProcedureReturn 2 ; nicht positiv definit
EndIf
mat(i,i) = SQRT(sum)
For j = i + 1 To n-1 Step 1
sum = mat(j,i)
For k = 0 To i-1 Step 1
sum - mat(i,k) * mat(j,k)
Next k
mat(j,i) = sum / mat(i,i)
Next j
Next i
ProcedureReturn 0
EndProcedure
;*====================================================================*
;* *
;* Num_ChoSol bestimmt die Loesung x des linearen Gleichungssystems *
;* B * x = b mit der unteren Dreieckmatrix B wie sie als Ausgabe *
;* von chodec bestimmt wird. *
;* *
;*====================================================================*
;* *
;* Eingabeparameter: *
;* ================ *
;* n Long n ( n > 0 ) *
;* Dimension von lmat, Anzahl der Komponenten *
;* des b-Vektors u. des Loesungsvektors x. *
;* lmat Double lmat(n,n) *
;* untere Dreieckmatrix, wie sie von chodec als Aus- *
;* gabe geliefert wird. *
;* b Double b(n) *
;* Rechte Seite des Gleichungssystems. *
;* *
;* Ausgabeparameter: *
;* ================ *
;* x Double x(n) *
;* Loesungsvektor des Systems. *
;* *
;* Rueckgabewert: *
;* ============= *
;* = 0 alles ok *
;* = 1 Zerlegungsmatrix unzulaessig oder n < 1 *
;* *
;*====================================================================*
Procedure.l Num_ChoSol (n.l, Array lmat.d(2), Array b.d(1), Array x.d(1))
Protected j.l, k.l
Protected sum.d
If n < 1
ProcedureReturn 1
EndIf
If lmat(0,0) = 0.0
ProcedureReturn 1 ; Unzulaessige Zerlegungsmatrix
EndIf
x(0) = b(0) / lmat(0,0) ; Vorwaertselimination
For k = 1 To n-1 Step 1
sum = 0.0
For j = 0 To k-1 Step 1
sum + lmat(k,j) * x(j)
Next j
If lmat(k,k) = 0.0
ProcedureReturn 1
EndIf
x(k) = (b(k) - sum) / lmat(k,k)
Next k
x(n-1) / lmat(n-1,n-1) ; Rueckwaertselimination
For k = n - 2 To 0 Step -1
sum = 0.0
For j = k + 1 To n-1 Step 1
sum + lmat(j,k) * x(j)
Next j
x(k) = (x(k) - sum) / lmat(k,k)
Next k
ProcedureReturn 0
EndProcedure
;*====================================================================*
;* *
;* Die Funktion cholesky dient zur Loesung eines linearen *
;* Gleichungssystems: mat * x = b *
;* mit der n x n Koeffizientenmatrix mat und der rechten Seite b *
;* nach dem Cholesky-Verfahren. *
;* *
;* Die Eingabematrix muss symmetrisch und positiv definit sein, *
;* d.h. fuer alle n Vektoren y != 0 muss gelten: *
;* T *
;* y * mat * y > 0. *
;* *
;* cholesky arbeitet nur auf der unteren Dreieckmatrix incl. der *
;* Diagonalen, so dass gegenueber dem Gaussverfahren eine erhebliche *
;* Reduktion des Speicherplatzes erreicht wird; es genuegt daher *
;* diesen Teil der Matrix an cholesky zu uebergeben. *
;* *
;*====================================================================*
;* *
;* Anwendung: *
;* ========= *
;* *
;* Lineare Gleichungssysteme mit n x n Koeffizientenmatrix, *
;* die symmetrisch und positiv definit ist. *
;* *
;*====================================================================*
;* *
;* Steuerparameter: *
;* =============== *
;* mod Long mod; *
;* Aufrufart von cholesky: *
;* = 0 Bestimmung der Zerlegungsmatrix und Berechnung *
;* der Loesung des Gleichungssystems. *
;* = 1 Nur Berechnung der Zerlegungsmatrix; wird auf *
;* mat ueberspeichert. *
;* = 2 Nur Loesung des Gleichungssystems; zuvor muss je- *
;* doch die Zerlegungsmatrix bestimmt sein. Diese *
;* Aufrufart wird verwendet, falls bei gleicher *
;* Matrix lediglich die rechte Seite des Systems vari- *
;* iert, z. B. zur Berechnung der Inversen. *
;* *
;* Eingabeparameter: *
;* ================ *
;* n Long n; ( n > 0 ) *
;* Dimension von mat, Anzahl der Komponenten *
;* des b-Vektors u. des Loesungsvektors x. *
;* mat Double mat(n,n) *
;* mod = 0, 1: Matrix des Gleichungssystems. *
;* mod = 2 : Zerlegungsmatrix. *
;* b Double b(n) ( bei mod = 0, 2 ) *
;* Rechte Seite des Gleichungssystems. *
;* *
;* Ausgabeparameter: *
;* ================ *
;* mat Double mat(n,n) ( bei mod = 0, 1 ) *
;* Dekompositionsmatrix, die die Zerlegung von *
;* mat enthaelt. *
;* x Double x(n); ( bei mod = 0, 2 ) *
;* Loesungsvektor des Systems. *
;* *
;* Rueckgabewert: *
;* ============= *
;* = 0 alles ok *
;* = 1 n < 1 gewaehlt oder falsche Eingabeparameter *
;* = 2 Matrix ist nicht positiv definit *
;* = 3 Falsche Aufrufart *
;* *
;*====================================================================*
;* *
;* Benutzte Funktionen: *
;* =================== *
;* Num_ChoDec() : Bestimmt die Dekomposition *
;* Num_ChoSol() : Loest das lineare Gleichungssystem *
;* *
;*====================================================================*
Procedure Num_Cholesky (mod.l,n.l,Array mat.d(2), Array b.d(1), Array x.d(1))
Protected rc.l
If n < 1
ProcedureReturn 1
EndIf
If mod = 0 ; Zerlegung bestimmen und Gleichungssystem loesen
rc = Num_ChoDec (n, mat())
If rc = 0
ProcedureReturn Num_ChoSol (n, mat(), b(), x())
Else
ProcedureReturn rc
EndIf
ElseIf mod = 1 ; Nur Zerlegung berechnen
ProcedureReturn Num_ChoDec (n, mat())
ElseIf mod = 2 ; Nur Loesung bestimmen
ProcedureReturn Num_ChoSol (n, mat(), b(), x())
EndIf
ProcedureReturn 3 ; Falsche Aufrufart
EndProcedure
;************************************************************************
;* die Koeffizienten eines Ausgleichspolynoms n-ten Grades nach der *
;* diskreten Fehlerquadratmethode von Gauss berechnen. *
;* Das Verfahren beruht auf den Gaussschen Normalgleichungen, die mit *
;* dem Choleskyverfahren geloest werden: Die Koeffizienten c[j] des *
;* Vektors *
;* *
;* ( w[0] * (c[0] + c[1] * x[0]^1 + ... + c[n] * x[0]^n) ) *
;* ( w[1] * (c[0] + c[1] * x[1]^1 + ... + c[n] * x[1]^n) ) *
;* ( ... ... ) *
;* ( w[m] * (c[0] + c[1] * x[m]^1 + ... + c[n] * x[m]^n) ) *
;* *
;* sollen so berechnet werden, dass er der rechten Seite *
;* *
;* ( w[0] * y[0] ) *
;* ( w[1] * y[1] ) *
;* ( ... ) *
;* ( w[m] * y[m] ) *
;* *
;* in der Euklidischen Norm moeglichst nahe kommt. *
;* c ergibt sich dann als Loesung des linearen Gleichungssystems *
;* a * c = b (Normalgleichungen) *
;* mit *
;* a[i][j] = w[0] * x[0]^(j+i) + ... + w[m] * x[m]^(j+i), *
;* b[i] = w[0] * y[0] * x[0]^i + ... + w[m] * y[m] * x[m]^i *
;* fuer i,j=0(1)n. *
;* *
;* Eingabeparameter: *
;* ================= *
;* m: Nummer der letzten Stuetzstelle *
;* n: Hoechstgrad des algebraischen Ausgleichspolynoms *
;* x: [0..m]-Vektor mit den Stuetzstellen *
;* y: [0..m]-Vektor mit den Funktionswerten zu den Stuetzstellen *
;* w: [0..m]-Vektor mit den Gewichten *
;* *
;* Ausgabeparameter: *
;* ================= *
;* c: [0..n]-Vektor mit den Koeffizienten des Ausgleichspolynoms *
;* *
;* Funktionswert: *
;* ============== *
;* 0: kein Fehler *
;* 1: Fehler in den Eingabedaten: n zu klein oder m zu klein *
;* 2: Fehler in choly() *
;* *
;************************************************************************
Procedure Num_ApproxPoly (n.l,m.l,Array x.d(1),Array y.d(1), Array w.d(1),Array c.d(1))
Dim a.d(n,n) ; die [0..n,0..n]-Matrix des Gleichungssystems
Dim b.d(n) ; [0..n]-Vektor mit der rechten Seite des Gleichungssystems
Protected summand.d ; Hilfsvariable zur Berechnung eines Matrixelements
Protected i.l, j.l ; Laufvariablen
If n < 0 Or m < n ; n oder m zu klein?
ProcedureReturn 1
EndIf
; die erste Spalte, die letzte Zeile der Matrix und die rechte Seite des Gleichungssystems berechnen
For i = 0 To n Step 1
a(i,0) = 0.0
a(n,i) = 0.0
b(i) = 0.0
Next i
For j = 0 To m Step 1
summand = w(j)
For i = 0 To n-1 Step 1
a(i,0) + summand
b(i) + summand * y(j)
summand * x(j)
Next i
a(n,0) + summand
b(n) + summand * y(j)
For i = 1 To n Step 1
summand * x(j)
a(n,i) + summand
Next i
Next j
;die Matrix vervollstaendigen, indem immer wieder Zeile i + 1 versetzt in Zeile i kopiert wird
For i = n - 1 To 0 Step -1
For j = 0 To n-1 Step 1
a(i,j + 1) = a(i + 1,j)
Next j
Next i
i = Num_Cholesky(0, n + 1, a(), b(), c()); ; das Gleichungssystem loesen
If i <> 0 ; Fehler beim Loesen des Gleichungsystems?
ProcedureReturn 2
EndIf
ProcedureReturn 0
EndProcedure
;=======================================================================================
;- Demo
#Image1 = 1
#Window = 2
#ImageGad = 3
#Button = 4
#Option1 = 5
#Option2 = 6
#MaxX = 800
#MaxY = 700
#MaxAnz = 1000
Dim X.d(#MaxAnz)
Dim Y.d(#MaxAnz)
Dim W.d(#MaxAnz)
Dim C.d(#MaxAnz)
Define i.l
Define Anzahl.l = 0
Define Grad.l = 0
Define Event.l
Define Aktion.l = #Option1
Define Input$
Define PosX.l
Define PosY.l
Define Xmin.d
Define Xmax.d
Define XA.l
Define XE.l
Define YA.l
Define YE.l
CreateImage(#Image1,#MaxX,#MaxY)
If OpenWindow(#Window,0,0,#MaxX+10,#MaxY+55,"Polynomial regression",#PB_Window_ScreenCentered)
ImageGadget (#ImageGad, 5,10, #MaxX,#MaxY,ImageID(#Image1))
ButtonGadget(#Button,5,#MaxY+20,100,25, "Quit")
OptionGadget(#Option1,110,#MaxY+15,100,20,"Add Points")
OptionGadget(#Option2,110,#MaxY+35,100,20,"Polynomial regression")
StartDrawing(ImageOutput(#Image1))
Box(0,0,#MaxX,#MaxY,RGB(50,50,50))
StopDrawing()
SetGadgetState(#ImageGad,ImageID(#Image1))
SetGadgetState(#Option1,1)
Repeat
Event = WaitWindowEvent()
Select Event
Case #PB_Event_Gadget
Select EventGadget()
Case #ImageGad
If Aktion = #Option1 And EventType() = #PB_EventType_LeftClick And Anzahl < #MaxAnz
PosX = WindowMouseX(#Window) - GadgetX(#ImageGad)
PosY = WindowMouseY(#Window) - GadgetY(#ImageGad)
X(Anzahl) = PosX
Y(Anzahl) = PosY
W(Anzahl) = 1.0
Anzahl + 1
StartDrawing(ImageOutput(#Image1))
Box(0,0,#MaxX,#MaxY,RGB(50,50,50))
For i = 0 To Anzahl-1
Box (Int(X(i)) - 2,Int(Y(i)) - 2,4,4,RGB(250,250,0))
Next
StopDrawing()
SetGadgetState(#ImageGad,ImageID(#Image1))
EndIf
Case #Button
Event = #PB_Event_CloseWindow
Case #Option1
If Aktion <> #Option1
Anzahl = 0
Grad = 0
StartDrawing(ImageOutput(#Image1))
Box(0,0,#MaxX,#MaxY,RGB(50,50,50))
StopDrawing()
SetGadgetState(#ImageGad,ImageID(#Image1))
Aktion = #Option1
EndIf
Case #Option2
If Aktion <> #Option2
If Anzahl >= 2
Repeat
Grad = Anzahl-1
Input$ = InputRequester("Degree of polynomial", "(min.:1/max."+Str(Grad)+"):", Str(Grad))
Grad = Val(Input$)
Until Grad>0 And Grad < Anzahl
i = Num_ApproxPoly(Grad,Anzahl-1,X(),Y(),W(),C())
If i > 0
MessageRequester("Error!","Error (" + Str(i) + ") during the calculation!")
StartDrawing(ImageOutput(#Image1))
Box(0,0,#MaxX,#MaxY,RGB(50,50,50))
StopDrawing()
SetGadgetState(#ImageGad,ImageID(#Image1))
Aktion = #Option1
Anzahl = 0
SetGadgetState(#Option1,1)
Else
StartDrawing(ImageOutput(#Image1))
Box(0,0,#MaxX,#MaxY,RGB(50,50,50))
Xmin = #MaxX
Xmax = 0.0
For i = 0 To Anzahl-1
Box (Int(X(i)) - 2,Int(Y(i)) - 2,4,4,RGB(250,250,0))
If X(i) < Xmin
Xmin = X(i)
EndIf
If X(i) > Xmax
Xmax = X(i)
EndIf
Next
XA = Int(Xmin)
XE = Int(Xmax)
For i=XA To XE Step 1
YE = Num_Horner(Grad,C(),i)
If i>XA
LineXY(i-1,YA,i,YE,RGB(255,255,255))
EndIf
YA = YE
Next i
StopDrawing()
SetGadgetState(#ImageGad,ImageID(#Image1))
Aktion = #Option2
EndIf
Else
MessageRequester("Error!","To less points")
SetGadgetState(#Option1,1)
EndIf
EndIf
EndSelect
EndSelect
Until Event = #PB_Event_CloseWindow
EndIf
CloseWindow(#Window)
PB 6.01 ― Win 10, 21H2 ― Ryzen 9 3900X, 32 GB ― NVIDIA GeForce RTX 3080 ― Vivaldi 6.0 ― www.unionbytes.de
Lizard - Script language for symbolic calculations and more ― Typeface - Sprite-based font include/module
Lizard - Script language for symbolic calculations and more ― Typeface - Sprite-based font include/module
- zxretrosoft
- Enthusiast
- Posts: 169
- Joined: Wed May 15, 2013 8:26 am
- Location: Czech Republic, Prague
- Contact:
Re: Cubic regression
Yes, this is an eventuality.
But it is necessary to solve it numerically.
x=1 y=5
x=2 y=7
x=3 y=11
x=4 y=20
Result (by the calculator):
A=0.5
B=-2
C=4.5
D=1.999999999
How did the calculator calculate it?
There is a formula for ABCD.
But it is necessary to solve it numerically.
x=1 y=5
x=2 y=7
x=3 y=11
x=4 y=20
Result (by the calculator):
A=0.5
B=-2
C=4.5
D=1.999999999
How did the calculator calculate it?
There is a formula for ABCD.
- NicTheQuick
- Addict
- Posts: 1227
- Joined: Sun Jun 22, 2003 7:43 pm
- Location: Germany, Saarbrücken
- Contact:
Re: Cubic regression
Was this a question? Then no. You only have to solve it numerically if you have less input values than output values.zxretrosoft wrote:But it is necessary to solve it numerically.
The english grammar is freeware, you can use it freely - But it's not Open Source, i.e. you can not change it or publish it in altered way.
- zxretrosoft
- Enthusiast
- Posts: 169
- Joined: Wed May 15, 2013 8:26 am
- Location: Czech Republic, Prague
- Contact:
Re: Cubic regression
Of course.
The answer is very very simple.
The quadratic regression formulas for ABC are:
a = { [ Σ x2 y * Σ xx ] - [Σ xy * Σ xx2 ] } / { [ Σ xx * Σ x2x 2] - [Σ xx2 ]2 }
b = { [ Σ xy * Σ x2x2 ] - [Σ x2y * Σ xx2 ] } / { [ Σ xx * Σ x2x 2] - [Σ xx2 ]2 }
c = [ Σ y / n ] - { b * [ Σ x / n ] } - { a * [ Σ x 2 / n ] }
The cubic regression formulas for ABCD are:
a =
b =
c =
d =
The answer is very very simple.
The quadratic regression formulas for ABC are:
a = { [ Σ x2 y * Σ xx ] - [Σ xy * Σ xx2 ] } / { [ Σ xx * Σ x2x 2] - [Σ xx2 ]2 }
b = { [ Σ xy * Σ x2x2 ] - [Σ x2y * Σ xx2 ] } / { [ Σ xx * Σ x2x 2] - [Σ xx2 ]2 }
c = [ Σ y / n ] - { b * [ Σ x / n ] } - { a * [ Σ x 2 / n ] }
The cubic regression formulas for ABCD are:
a =
b =
c =
d =
Re: Cubic regression
Have you read the article?
Have you seen my code?
Of course you can use a very very long formula for A,B,C,D or you can just solve/calculate it with the matrix algorithm. The result of A, B, C and D is in the Array c.d()
Your formula for degree 2 or 3 is nothing else like this matrix formula, just rewritten.
Btw:
Have you seen my code?
Of course you can use a very very long formula for A,B,C,D or you can just solve/calculate it with the matrix algorithm. The result of A, B, C and D is in the Array c.d()
Your formula for degree 2 or 3 is nothing else like this matrix formula, just rewritten.
Which calculator?How did the calculator calculate it?
Btw:
Code: Select all
x(0) = 1 : y(0) = 5
x(1) = 2 : y(1) = 7
x(2) = 3 : y(2) = 11
x(3) = 4 : y(3) = 20
w(0) = 1
w(1) = 1
w(2) = 1
w(3) = 1
Num_ApproxPoly (3, 3, x(), y(), w(), c())
Debug "D = "+c(0)
Debug "C = "+c(1)
Debug "B = "+c(2)
Debug "A = "+c(3)
D = 2
C = 4.5
B = -2
A = 0.5
Last edited by STARGÅTE on Wed Jul 17, 2019 12:51 pm, edited 2 times in total.
PB 6.01 ― Win 10, 21H2 ― Ryzen 9 3900X, 32 GB ― NVIDIA GeForce RTX 3080 ― Vivaldi 6.0 ― www.unionbytes.de
Lizard - Script language for symbolic calculations and more ― Typeface - Sprite-based font include/module
Lizard - Script language for symbolic calculations and more ― Typeface - Sprite-based font include/module
- zxretrosoft
- Enthusiast
- Posts: 169
- Joined: Wed May 15, 2013 8:26 am
- Location: Czech Republic, Prague
- Contact:
Re: Cubic regression
Yes, thank you, I'm studying the code It's difficult for me because I can't speak English. It takes me to understand something.STARGÅTE wrote:Have you read the article?
Have you seen my code?
Of course you can use a very very long formula for A,B,C,D or you can just solve/calculate it with the matrix algorithm. The result of A, B, C and D is in the Array c.d()
Your formula for degree 2 or 3 is nothing else like this matrix formula, just rewritten.
Which calculator?How did the calculator calculate it?
In the code, where can I enter those XY values?
P.S. Calculator Casio fx9860.
Re: Cubic regression
Copy the code (without the DEMO) and add this:
Code: Select all
; Enter here all points x1,y1 ; x1,y2 ; ... :
Define PointList.s = "1,5 ; 2,7 ; 3,11 ; 4,20"
Define Points.i = CountString(PointList, ";")
Define Dim X.d(Points)
Define Dim Y.d(Points)
Define Dim W.d(Points)
Define Dim C.d(3) ; Degree of polynom
Define I
For I = 0 To Points
X(I) = ValD(StringField(StringField(PointList, I+1, ";"), 1, ","))
Y(I) = ValD(StringField(StringField(PointList, I+1, ";"), 2, ","))
W(I) = 1.0
Next
If Num_ApproxPoly(Points, 3, X(), Y(), W(), C()) = 0
Debug "A = "+C(3)
Debug "B = "+C(2)
Debug "C = "+C(1)
Debug "D = "+C(0)
Else
Debug "Error"
EndIf
A = 0.5
B = -2
C = 4.5
D = 2
PB 6.01 ― Win 10, 21H2 ― Ryzen 9 3900X, 32 GB ― NVIDIA GeForce RTX 3080 ― Vivaldi 6.0 ― www.unionbytes.de
Lizard - Script language for symbolic calculations and more ― Typeface - Sprite-based font include/module
Lizard - Script language for symbolic calculations and more ― Typeface - Sprite-based font include/module
- zxretrosoft
- Enthusiast
- Posts: 169
- Joined: Wed May 15, 2013 8:26 am
- Location: Czech Republic, Prague
- Contact:
Re: Cubic regression
Yes it works right! Your code is great and it's the solution I'm looking for. I'll still use the code. Thank you so much STARGÅTE!
But I still search for 4 ABCD equations that I can't find. I don't care if they're long, but I'd like to know them.
But I still search for 4 ABCD equations that I can't find. I don't care if they're long, but I'd like to know them.
- NicTheQuick
- Addict
- Posts: 1227
- Joined: Sun Jun 22, 2003 7:43 pm
- Location: Germany, Saarbrücken
- Contact:
Re: Cubic regression
I don't think this will help you in any way. Just understand that you better should use matrix calculations with more than 3 values.zxretrosoft wrote:But I still search for 4 ABCD equations that I can't find. I don't care if they're long, but I'd like to know them.
The english grammar is freeware, you can use it freely - But it's not Open Source, i.e. you can not change it or publish it in altered way.
Re: Cubic regression
using Maple the symbolic coefficients are
you could probably rewrite them using sums and therefore making the expressions shorter and more aesthetically pleasing, I leave that up to you
one could achieve the same result using the free Maxima CAS https://wxmaxima-developers.github.io/wxmaxima/
maxima code follows
Code: Select all
A = y1/(x1^3-x1^2*x2-x1^2*x3-x1^2*x4+x1*x2*x3+x1*x2*x4+x1*x3*x4-x2*x3*x4)-y2/(x1*x2^2-x1*x2*x3-x1*x2*x4+x1*x3*x4-x2^3+x2^2*x3+x2^2*x4-x2*x3*x4)+y3/(x1*x2*x3-x1*x2*x4-x1*x3^2+x1*x3*x4-x2*x3^2+x2*x3*x4+x3^3-x3^2*x4)-y4/(x1*x2*x3-x1*x2*x4-x1*x3*x4+x1*x4^2-x2*x3*x4+x2*x4^2+x3*x4^2-x4^3)
B = -(x2+x3+x4)*y1/(x1^3-x1^2*x2-x1^2*x3-x1^2*x4+x1*x2*x3+x1*x2*x4+x1*x3*x4-x2*x3*x4)+(x1+x3+x4)*y2/(x1*x2^2-x1*x2*x3-x1*x2*x4+x1*x3*x4-x2^3+x2^2*x3+x2^2*x4-x2*x3*x4)-(x1+x2+x4)*y3/(x1*x2*x3-x1*x2*x4-x1*x3^2+x1*x3*x4-x2*x3^2+x2*x3*x4+x3^3-x3^2*x4)+(x1+x2+x3)*y4/(x1*x2*x3-x1*x2*x4-x1*x3*x4+x1*x4^2-x2*x3*x4+x2*x4^2+x3*x4^2-x4^3)
C = (x2*x3+x2*x4+x3*x4)*y1/(x1^3-x1^2*x2-x1^2*x3-x1^2*x4+x1*x2*x3+x1*x2*x4+x1*x3*x4-x2*x3*x4)-(x1*x3+x1*x4+x3*x4)*y2/(x1*x2^2-x1*x2*x3-x1*x2*x4+x1*x3*x4-x2^3+x2^2*x3+x2^2*x4-x2*x3*x4)+(x1*x2+x1*x4+x2*x4)*y3/(x1*x2*x3-x1*x2*x4-x1*x3^2+x1*x3*x4-x2*x3^2+x2*x3*x4+x3^3-x3^2*x4)-(x1*x2+x1*x3+x2*x3)*y4/(x1*x2*x3-x1*x2*x4-x1*x3*x4+x1*x4^2-x2*x3*x4+x2*x4^2+x3*x4^2-x4^3)
D = -x2*x3*x4*y1/(x1^3-x1^2*x2-x1^2*x3-x1^2*x4+x1*x2*x3+x1*x2*x4+x1*x3*x4-x2*x3*x4)+x4*x3*x1*y2/(x1*x2^2-x1*x2*x3-x1*x2*x4+x1*x3*x4-x2^3+x2^2*x3+x2^2*x4-x2*x3*x4)-x1*x2*x4*y3/(x1*x2*x3-x1*x2*x4-x1*x3^2+x1*x3*x4-x2*x3^2+x2*x3*x4+x3^3-x3^2*x4)+x1*x2*x3*y4/(x1*x2*x3-x1*x2*x4-x1*x3*x4+x1*x4^2-x2*x3*x4+x2*x4^2+x3*x4^2-x4^3)
one could achieve the same result using the free Maxima CAS https://wxmaxima-developers.github.io/wxmaxima/
maxima code follows
Code: Select all
load ("eigen")$
m:matrix([1,x1,x1^2,x1^3],[1,x2,x2^2,x2^3], [1,x3,x3^2,x3^3], [1,x4,x4^2,x4^3]);
y:matrix([y1],[y2], [y3], [y4]);
b:invert_by_lu(m).y$
A:ratsimp(b[4][1]);
B:ratsimp(b[3][1]);
C:ratsimp(b[2][1]);
D:ratsimp(b[1][1]);
subst([x1=1, x2=2,x3=3,x4=4,y1=5,y2=7,y3=11,y4=20],A);
subst([x1=1, x2=2,x3=3,x4=4,y1=5,y2=7,y3=11,y4=20],B);
subst([x1=1, x2=2,x3=3,x4=4,y1=5,y2=7,y3=11,y4=20],C);
subst([x1=1, x2=2,x3=3,x4=4,y1=5,y2=7,y3=11,y4=20],D);
- zxretrosoft
- Enthusiast
- Posts: 169
- Joined: Wed May 15, 2013 8:26 am
- Location: Czech Republic, Prague
- Contact:
Re: Cubic regression
Yes jack! This is exactly it!
But it seems to me that it's just for 4 input values. Can we generalize it anyway? (i.e. for any number of values) ?
But it seems to me that it's just for 4 input values. Can we generalize it anyway? (i.e. for any number of values) ?
- NicTheQuick
- Addict
- Posts: 1227
- Joined: Sun Jun 22, 2003 7:43 pm
- Location: Germany, Saarbrücken
- Contact:
Re: Cubic regression
The generalization is the matrix calculation STARGÅTE posted.
The english grammar is freeware, you can use it freely - But it's not Open Source, i.e. you can not change it or publish it in altered way.