Algoritmo per trovare per tentativi, con R, la funzione di un certo tipo che "meglio" approssima certi dati  (viene calcolata la somma dei quadrati degli scarti tra i dati (input,output) e le immagini degli input mediante la funzione).
  Occorre scegliere un intervallo iniziale ampio [m1,M1] per il parametro incognito U ( o 2 o 3 o 4 intervalli inziali ampi [m1,M1], [m2,M2], [m3,M3], [m4,M4] per i parametri incogniti U, V, W, Z) che contengano sicuramente le soluzioni (occorre sceglierli tenendo conto dei punti dati e del tipo di funzione con cui li si vogliono approssimare), e dopo il comando F_RUN1 (o F_RUN2, F_RUN3, F_RUN4).  Poi si possono man mano restringere gli intervalli.
  Occorre scegliere un numero di passi tanto più alto quanti di più sono i parametri.
  Ovviamente l'algoritmo deve essere modificato a seconda del tipo di funzione che si vuole cercare.
  Se proseguissi gli studi nell'ambito della matematica applicata vedresti che esistono delle tecniche abbastanza sofisticate per trovare le curve per via teorica.
  Per usare il programma prima caricare source("http://macosa.dima.unige.it/r.R") e poi le righe seguenti, dopo aver manipolato quelle in blu. Vedi gli esempi  uno,..., sei.

x = c( ... ); y = c( ... )
# Scegli e definisci il tipo di funzione ...
yy1 = function(U,x) ...; Fun = function(x) ... # con U0 al posto di U
# o:
yy2 = function(U,V,x) ...; Fun = function(x) ...
# o:
yy3 = function(U,V,W,x) ...; Fun = function(x) ...
# o:
yy4 = function(U,V,W,Z,x) ...; Fun = function(x) ...
# In U0 (V0,W0,Z0) vengono messi da  F_RUN1 (… F_RUN4)  i valori che minimizzano lo scarto
# misurato come somma dei quadrati delle differenze tra le ordinate e i corrispondenti
# valori di Fun(x), e in  dif  lo scarto ottenuto, di cui stampo la radice quadrata.
# All'inizio scelgo uno scarto grande, che verrà minimizzato.
dif = 1e300
# Inzio scegliendo intervalli in cui stanno di sicuro i parametri incogniti:
m1=…; M1=…; F_RUN1(x,y,1e5); cat(U0," d:",NUMBER(dif^0.5),"\n")
# o:
# ...
m1=…;M1=…; m2=…;M2=…; m3=…;M3=…; m4=…;M4=…; F_RUN4(x,y,1e5); cat(U0,V0,W0,Z0," d:",NUMBER(dif^0.5),"\n")
# Proseguo restringendo gli intervalli; eventualmente aumento il n. delle prove
# ...
# Quando voglio traccio il grafico della funzione man mano ottenuta
# Alla fine traccio il grafico della funzione su cui ci si è stabilizzati

Esempi

Uno
Voglio approssimare la relazione tra due variabili casuali G1 e G2 con una relazione del tipo G2 = k G1 + h.
Devo determinare k ed h. In tre esperimenti ottengo le coppie (1.6,18), (3.6,31), (4.8,48).

G1 = c(1.6, 3.6, 4.8)
G2 = c(18,  31,  48)
range(G1); range(G2)
# 1.6  4.8    18  48
Plane(0,5, 0,50); POINT(G1,G2,"brown")
yy2 = function(U,V,x) U+V*x; Fun = function(x) U0+V0*x
dif = 1e300
m1=-10; M1=10; m2=0; M2=50; F_RUN2(G1,G2,1e5); cat(U0,V0, " d:",NUMBER(dif^0.5),"\n")
#2.06595 9.074395  d: 4.64685263153777
# facciamo il grafico per vedere che cosa abbiamo ottenuto
graph1(Fun,0,5, "red")
# In questo caso (lineare) proseguendo ulteriormente non miglioro in modo significativo:
m1=1; M1=3; m2=8; M2=10; F_RUN2(G1,G2,1e5); cat(U0,V0, " d:",NUMBER(dif^0.5),"\n")
#2.068302 9.079689  d: 4.64670394481248
m1=2; M1=2.08; m2=9; M2=9.1; F_RUN2(G1,G2,1e5); cat(U0,V0, " d:",NUMBER(dif^0.5),"\n")
#2.061353 9.081599  d: 4.64670170566259    2.061 + 9.082*x
# volendo proseguire:
m1=2.06; M1=2.065; m2=9.075; M2=9.085; F_RUN2(G1,G2,1e5); cat(U0,V0, " d:",NUMBER(dif^0.5),"\n")
#2.061236 9.081633  d: 4.64670170498558
m1=2.061; M1=2.062; m2=9.081; M2=9.082; F_RUN2(G1,G2,1e5); cat(U0,V0, " d:",NUMBER(dif^0.5),"\n")
#2.061222 9.081633  d: 4.6467017049404     2.0612 + 9.0816*x
graph1(Fun,0,5, "brown")
                
# Ovviamente in questo caso avrei potuto usare la tecnica di regressione:
regression1(G1,G2)
# 9.082 * x + 2.061
Due
So che una grandezza y legata ad una grandezza x da una relazione del tipo y = u*exp(v*x)+w. Conosco gli esiti sperimentali (1.186, 0.955, 0.786, 0.648, 0.568, 0.500), a cui non so assegnare una precisione, associati ai valori 1, 2, 3, 4, 5, 6.  Devo determinare u, v e w.

x=c(1,2,3,4,5,6); y=c(1.186,0.955,0.786,0.648,0.568,0.500)
Plane(0,7, 0,1.3); POINT(x,y, "brown")
yy3 = function(U,V,W,x) U*exp(V*x)+W; Fun = function(x) U0*exp(V0*x)+W0
dif = 1e300
m1=0; M1=10; m2=-100; M2=0; m3=0; M3=1/2; F_RUN3(x,y,1e5); cat(U0,V0,W0," d:",NUMBER(dif^0.5),"\n")
1.091296 -0.3467326 0.3593317  d: 0.0881476985674848 
graph1(Fun,0,7, "red")
m1=0.5; M1=1.5; m2=-1; M2=0; m3=0.2; M3=0.4; F_RUN3(x,y,1e5); cat(U0,V0,W0," d:",NUMBER(dif^0.5),"\n")
1.19855 -0.3138632 0.3137454  d: 0.0113071865924665
...
m1=1.19; M1=1.2; m2=-0.35; M2=-0.3; m3=0.3; M3=0.35; F_RUN3(x,y,1e5); cat(U0,V0,W0," d:",NUMBER(dif^0.5),"\n")
1.190008 -0.3145388 0.318526  d: 0.0104177052956292 
graph1(Fun,0,7, "brown")
# y = 1.190*exp(-0.3145*x) + 0.3185

Tre
Un peso oscilla su una molla verticale. So che la sua altezza h (in cm) al variare del tempo t (in decimi di secondo) è h = U·cos(ωt)+V·sin(ωt). Riesco a misurare ω e trovo che vale 1.03 s−1. Attraverso una ripresa filmata trovo che a t pari, in ordine, a 1, 3, 5, 7, 9 corrispondono, approssimativamente, i valori di h 3, −16, 6, 9, −8. Stima quanto valgono (in cm) U e V.

t = c(1,   3, 5, 7,  9)
h = c(3, -16, 6, 9, -8)
range(t);range(h)
#  1 9   -16   9
BF=5; HF=3
Plane(0,10, -20,15); POINT(t,h, "brown")
omega = 1.03; yy2 = function(U,V,x) U*cos(omega*x)+V*sin(omega*x)
Fun = function(x) U0*cos(omega*x)+V0*sin(omega*x)
dif=1e300
m1=-100; M1=100; m2=-100; M2=100; F_RUN2(t,h,1e5); cat(U0,V0," d:",NUMBER(dif^0.5),"\n")
# 11.89429 -0.4119553  d: 6.64664848751153
m1=-50; M1=50; m2=-10; M2=10; F_RUN2(t,h,1e5); cat(U0,V0," d:",NUMBER(dif^0.5),"\n")
# 12.15634 -0.7992009  d: 6.61732609665133
m1=10; M1=15; m2=-1; M2=-0.5; F_RUN2(t,h,1e5); cat(U0,V0," d:",NUMBER(dif^0.5),"\n")
# 12.10602 -0.7896324  d: 6.61680941853707
m1=12; M1=13; m2=-0.8; M2=-0.75; F_RUN2(t,h,1e5); cat(U0,V0," d:",NUMBER(dif^0.5),"\n")
# 12.10647 -0.7943228  d: 6.61680536155479
m1=12.1; M1=12.11; m2=-0.8; M2=-0.79; F_RUN2(t,h,1e5); cat(U0,V0," d:",NUMBER(dif^0.5),"\n")
# 12.10616 -0.7946137  d: 6.61680532529941
m1=12.106; M1=12.107; m2=-0.795; M2=-0.794; F_RUN2(t,h,1e5); cat(U0,V0," d:",NUMBER(dif^0.5),"\n")
# 12.10615 -0.794591  d: 6.61680532519233
           
# h = 12.11*cos(omega*t)-0.795*sin(omega*t)

Quattro
Qual è la funzione x → A·x^3+B·x^2+C·x il cui grafico meglio approssima i punti (−5,4), (−2,−1), (−1,−3), (1,5), (3,0), (4,−6)? Se non imponessimo che il grafico passi per (0,0) potremmo trovarla con regression3.

x = c(-5,-2,-1,1,3, 4); y = c(4, -1,-3,5,0,-6)
range(x) ;range(y)
# -5  4   -6  5
Plane(-6,5, -7,6); POINT(x,y, "brown")
yy3 = function(U,V,W,x) U*x^3+V*x^2+W*x; Fun = function(x) U0*x^3+V0*x^2+W0*x 
# Per non procedere a caso cerco di capire come deve essere fatto piu'
# o meno il polinomio (dai punti capisco che ha coefficiente direttivo
# negativo e parametro di grado 3 non troppo grande)
dif = 1e300
m1 = -10; M1 = 0; m2 = -1e3; M2 = 1e3; m3 = -1e3; M3 = 1e3; F_RUN3(x,y,2e5); cat(U0,V0,W0," d:",NUMBER(dif^0.5),"\n")
# -1.346109 -3.806864 12.93589  d: 100.149417394534 
graph1(Fun, -7,6, "red")
m1 = -10; M1 = 0; m2 = -1e3; M2 = 1e3; m3 = -1e3; M3 = 1e3; F_RUN3(x,y,2e5); cat(U0,V0,W0," d:",NUMBER(dif^0.5),"\n")
# -0.3605754 -1.075944 1.372042  d: 33.967719669027
graph1(Fun, -7,6, "seagreen")
m1 = -2; M1 = 0; m2 = -5; M2 = 0; m3 = 0; M3 = 10; F_RUN3(x,y,2e5); cat(U0,V0,W0," d:",NUMBER(dif^0.5),"\n")
# -0.1448519 -0.207785 1.950992  d: 4.62038929996649
# ...
m1 = -0.155; M1 = -0.152; m2 = -0.22; M2 = -0.19; m3 = 1.8; M3 = 2.1; F_RUN3(x,y,2e5); cat(U0,V0,W0," d:",NUMBER(dif^0.5),"\n")
# -0.1539919 -0.2125809 1.928769  d: 4.39882991977069 
graph1(Fun,-7,6, "brown")
# y = -0.1540*x^3 - 0.2126*x^2 + 1.929*x
                  

Cinque
Trovare la funzione del tipo x → A+B·x+C·sin(π·x/2) il cui grafico meglio approssima i punti (0,−1), (1,1), (2,1), (3,2), (4,3).

x = c(0,1,2,3,4); y = c(-1,1,1,2,3)
Plane(-1,5, -2,4); POINT(x,y, "brown")
yy3 = function(U,V,W,x) U+V*x+W*sin(pi/2*x); Fun = function(x) U0+V0*x+W0*sin(pi/2*x)
dif = 1e300
m1=-5; M1 = 0; m2 = -1e2; M2 = 1e2; m3 = -1e2; M3 = 1e2; F_RUN3(x,y,2e5); cat(U0,V0,W0," d:",NUMBER(dif^0.5),"\n")
# -0.6856656 0.9383885 0.8853998  d: 0.854635046012712
graph1(Fun, -2,6, "red")
m1=-2; M1 = 0; m2 = -10; M2 = 10; m3 = -10; M3 = 10; F_RUN3(x,y,2e5); cat(U0,V0,W0," d:",NUMBER(dif^0.5),"\n")
# -0.8583343 0.9894515 0.501509  d: 0.576825403584979 
m1=-1.5; M1 = 0; m2 = 0; M2 = 2; m3 = 0; M3 = 1; F_RUN3(x,y,2e5); cat(U0,V0,W0," d:",NUMBER(dif^0.5),"\n")
# -0.8161146 1.00398 0.4829206  d: 0.548950093971122
# ...
m1=-0.82; M1 = -0.78; m2 = 0.95; M2 = 1.05; m3 = 0.45; M3 = 0.55; F_RUN3(x,y,2e5); cat(U0,V0,W0," d:",NUMBER(dif^0.5),"\n")
# -0.7999081 1.00002 0.5004336  d: 0.547722952719513
graph1(Fun, -2,6, "brown")
# y = -0.7999+x+0.50004*sin(pi/2*x)  o  y = -0.80+x+sin(pi/2*x)/2
                  

Sei
Trovare la funzione del tipo x → A·xB il cui grafico meglio approssima i punti (0.9,1.6), (1.4,2.4), (2,3.3), (3.2,4.4), (4,5).

x = c(0.9,1.4,2,3.2,4); y = c(1.6,2.4,3.3,4.4,5)
Plane(0,5, 0,6)
POINT(x,y, "brown")
yy2 = function(U,V,x) U*x^V; Fun = function(x) U0*x^V0
dif = 1e300
m1 = 0; M1 = 10; m2 = 0; M2 = 10; F_RUN2(x,y,1e5); cat(U0,V0," d:",NUMBER(dif^0.5),"\n")
# 1.929549 0.7078059  d: 0.28739501222601 
m1 = 0; M1 = 5; m2 = 0; M2 = 5; F_RUN2(x,y,1e5); cat(U0,V0," d:",NUMBER(dif^0.5),"\n")
# 1.909829 0.7070296  d: 0.272840336883223 
m1 = 1.5; M1 = 2.5; m2 = 0.5; M2 = 1; F_RUN2(x,y,1e5); cat(U0,V0," d:",NUMBER(dif^0.5),"\n")
# 1.900816 0.711557  d: 0.272237173340978 
m1 = 1.7; M1 = 2.1; m2 = 0.6; M2 = 0.8; F_RUN2(x,y,1e5); cat(U0,V0," d:",NUMBER(dif^0.5),"\n")
# 1.898605 0.7124717  d: 0.27221748879405 
m1 = 1.85; M1 = 1.95; m2 = 0.66; M2 = 0.75; F_RUN2(x,y,1e5); cat(U0,V0," d:",NUMBER(dif^0.5),"\n")
# 1.898828 0.7123531  d: 0.272217042392632 
m1 = 1.88; M1 = 1.91; m2 = 0.69; M2 = 0.73; F_RUN2(x,y,1e5); cat(U0,V0," d:",NUMBER(dif^0.5),"\n")
# 1.898733 0.7123597  d: 0.272216960655014 
graph1(Fun, 0,5, "brown")
# y = 1.899*x^0.7124  o  y = 1.9*x^0.712