Interpolazione di funzioni vettoriali in Matlab

[Ianero]

Buonasera a tutti,
qualcuno sa se esiste una funzione che provvede ad effettuare una interpolazione ai minimi quadrati, di un modello non lineare che gli fornisco io, di una funzione che restituisce più di un valore d'uscita?
Non è la stessa cosa che fare un fit sulle singole funzioni che definiscono la funzione vettoriale.

Grazie in anticipo.

[DrCox]

interpolazione ai minimi quadrati

Non esiste. O interpoli (passando per tutti i punti) o fai una approssimazione/regressione (ad es con i minimi quadrati).
Quello che vuoi fare non puoi chiamarlo interpolazione.

Per interpolare, la cosa che secondo me si avvicina di piu' e' usare la funzione

Codice: Seleziona tutto

interpn

Per la regressione, prova a guardare qua: https://www.mathworks.com/help/stats/mvregress.html
oppure prova a farla tu a mano.
A me non vengono in mente altre funzioni gia' pronte

[Ianero]

Mi sono accorto in realtà di aver commesso un errore. Quello che devo fare è trovare la curva 3D che meglio approssima un set di punti, condizionando il problema alla seguente equazione parametrica per la curva da trovare:



dove i parametri da stimare sono i tre angoli , , , e inoltre , e .
Dunque direi, se non faccio errori, che io devo minimizzare l'errore quadratico medio su due curve separatamente, e , eliminando il parametro tempo:





Fin qui ti torna o vedi errori di impostazione del problema?

Se è giusto, allora in Matlab lo farei così:

Codice: Seleziona tutto

x_fittype = fittype('(v_0*sin(beta)+9.8/k)/(v_0*cos(beta))*( (x-d_0*cos(alpha))/(cos(phi)) ) + 9.8/(k^2)*log( (1-k*(x-d_0*cos(alpha))/(cos(phi)*v_0*cos(beta))) )','independent','x','dependent','z')
opt=fitoptions('Method','NonlinearLeastSquares');
%specifico in ordine: alpha,beta,d_0,k,phi,v_0
opt.Lower      = [-pi/4 pi/4-5*pi/180  2000 1e-4 pi-10*pi/180 100];
opt.StartPoint = [0 pi/4  3000 1e-3 pi 150];
opt.Upper      = [pi/4 pi/4+5*pi/180  4500 1e-2 pi+10*pi/180 200];
z_fit = fit(x_data,z_data,x_fittype,opt)


e analogamente poi anche per .
Dandogli in pasto una traiettoria con , la stima su avviene benissimo, mentre quella su restituisce:

Complex value computed by model function, fitting cannot continue.
Try using or tightening upper and lower bounds on coefficients.


che sicuramente dipende da quel nel logaritmo.
Come si potrebbe gestire questo problema? Lui sta ricevendo tanti valori di z e corrispondentemente tanti valori di .

Grazie.

[xyz]

Sulle regressioni non lineari hai provato a usare R (programma open source multi piattaforma) ? Ha nativo il calcolo dei minimi quadrati non lineari:

https://www.rdocumentation.org/packages ... topics/nls

non mi è mai servito per funzioni vettoriali ma per funzioni singole non lineari individua i coefficienti in modo efficiente, basta caricare i dati e impostare il modello con i valori iniziali dei coefficienti da individuare.

[Ianero]

Vado a dare un’occhiata, grazie.
Nel frattempo ho scoperto che per il mio scopo forse ciò che mi serve è il metodo dei minimi quadrati totali (TLS), devo approfondire anche questo.

[DrCox]

Quello che devo fare è trovare la curva 3D che meglio approssima un set di punti, condizionando il problema alla seguente equazione parametrica per la curva da trovare:



dove i parametri da stimare sono i tre angoli , , , e inoltre , e .

Hai un certo numero di punti nello spazio e vuoi trovare il valore dei parametri del tuo modello al fine di fittare tale nuvola di punti.
Quello che farei io è:
- definire una funzione costo da minimizzare: puoi considerare ad esempio la somma dei quadrati delle distanze, o altre metriche.
- usare una funzione di minimizzazione che ti restituisca i parametri
Per quanto riguarda il secondo punto, invece di usare

Codice: Seleziona tutto

fit

, userei

Codice: Seleziona tutto

fmincon

oppure, in caso di problemi più rognosi, spesso mi ha salvato usare

Codice: Seleziona tutto

simulannealbnd

[Ianero]

Provo, grazie.

[xyz]

...il metodo dei minimi quadrati totali (TLS)...

Sarebbe in inglese Total Least Squares. E' sempre un metodo per individuare una retta, le funzioni che hai come modello sono trigonometriche, difficilmente riesci a trovare una trasformazione biettiva adatta allo scopo.

Qui spiega una tecnica con R per una regressione non-lineare con un modello una funzione vettoriale:

https://stats.stackexchange.com/questio ... ssion-in-r

[marioursino]

Seguo. Ho un problema simile.

[DrCox]

Alla fine della fiera, identificare dei parametri altro non e' che un problema di ottimizzazione, per il quale devi trovare il minimo di una certa funzione ("funzione costo") da te definita.

Partiamo da un caso piu' semplice, di cui e' immediato visualizzare il risultato.
Consideriamo questa funzione:




Come fittare un insieme di punti usando questa funzione?

Ecco un esempio di codice (che verosimilmente senza modificarlo non funzionera' con i tuoi dati, e' giusto per darti l'idea di come approcciare il problema)

main.m

Codice: Seleziona tutto

%
% data..... %already defined
% t..... %already defined
%

%% Initialize variables

alpha=0;
beta=0;
phi=0;
k=0;
v=0;
d0=0;

%% Fitting

x0 = [alpha, beta, phi, k, v, d0];
options = optimset('MaxFunEvals',10000,'TolFun',1e-9,'TolX',1e-9);
[x] = fminsearch(@(x) cost_function(t, data, x), x0, options);

alpha=x(1);
beta=x(2);
phi=x(3);
k=x(4);
v=x(5);
d0=x(6);

%% check result

[ x,y ] = fitting_function( t, alpha, beta, phi, k, v, d0 );

figure
plot3(t,data(1,:),data(2,:),'-b','LineWidth',2)
hold on
plot3(t,x,y,'--r','LineWidth',2)
xlabel('t')
ylabel('x')
zlabel('y')
grid minor
legend('data','fit')


cost_function.m

Codice: Seleziona tutto

function [ cost ] = cost_function( t, data, x )

alpha=x(1);
beta=x(2);
phi=x(3);
k=x(4);
v=x(5);
d0=x(6);

[ x_var, y_var ] = fitting_function( t, alpha, beta, phi, k, v, d0 );

cost_x = sum( (x_var - data(1,:)).^2 );
cost_y = sum( (y_var - data(2,:)).^2 );

cost = mean([cost_x, cost_y]);

end


fitting_function.m

Codice: Seleziona tutto

function [ x,y ] = fitting_function( t, alpha, beta, phi, k, v, d0 )

x = v*cos(phi)*cos(beta)./k * (1-exp(-k*t))+d0*cos(alpha);
y = v*sin(phi)*cos(beta)./k * (1-exp(-k*t))+d0*sin(alpha);

end


Per il problema x,y,z,t che hai tu, devi solo aggiungere una dimensione.


Se il fitting non funziona, per problemi con dati piu' complessi, si puo':
- cambiare tecnica di ottimizzazione (come accennato sopra, ad esempio la simulated-annealing)
- cambiare gli initial guess
- definire una cost-function differente