MATLAB: How to determine constants using a linear solve of logarithmic transformations which minimize the sum of the residual squares for the below equations

constant determinationMATLABmultivariate linear regressionmultivariate power equationmultivariate regressionmvregress()regress

I have a multivariate power equation in the form y = a1 * (x1^a2) * (x2^a3) * (x3^a4) * (x4^a5) * (x5^a6) , where I have an array of values for each x-term and for the y-term that are each 19775×1. I need to determine the values for the constants a1, a2, a3, a4, and a5 and was told to use a linear solve of logarithmic transformations which minimize the sum of the residual squares. Is there a way to do this in MATLAB?

Attached is my code that I have thus far:

%% Non-dimensional empirical model
% fit used to create correlation:
% pi_theta = a1 .* (pi_y).^a2 .* (pi_z).^a3 .* (pi_twp).^a4 .* ...
%     (pi_tcb).^a5 .* (pi_alpha_cb).^a6;
% log10(pi_theta) = log10(a1) + a2 * log10(pi_y) + a3 * log10(pi_z)...
%     + a4 * log10(pi_twp) + a5 * log10(pi_tcb) + a6 * log10(pi_alpha_cb);
% solve constants using linear solve of logarithmic transformations of
% above equation.
% a1 = 
% a2 = 
% a3 = 
% a4 = 
% a5 = 
% a6 = 
X = [ones(size(pi_y)), log10(pi_y), log10(pi_z), log10(pi_twp), ...
    log10(pi_tcb), log10(pi_alpha_cb)];
Y = log10(pi_theta);
[b, bint,r,rint,stats] = mvregress(Y,X);

I keep getting this error, though:

Error using mvregress (line 456)
Covariance is not positive-definite.
Error in Dimensionless_parameters_analysis (line 266)
[b, bint,r,rint,stats] = mvregress(Y,X);

Best Answer

REPOSTED —

I decided to let ga have a go at your data and function.

This code :

D = load('X-Y.mat '); 
X = D.X ; 
Y = D.Y; 
objfcn = @(a,x) a(1) .* (x(:,1).^a(2)) .* (x(:,2).^a(3)) .* (x(:,3).^a(4)) .* (x(:,4).^a(5)) .* (x(:,5).^a(6 )); 
ftns = @(a) norm(Y - objfcn(a,X)); 
PopSz = 500 ; 
Parms = 6 ; 
opts = optimoptions('ga', 'PopulationSize',PopSz, 'InitialPopulationMatrix',randi(1E+4,PopSz,Parms)*1E-1, 'MaxGenerations',2E3, 'PlotFcn',@gaplotbestf, 'PlotInterval',1 ); 
[theta,fval,exitflag,output] = ga(ftns, Parms, [],[],[],[],zeros(Parms,1),Inf(Parms,1),[],[],opts )

converged fairly quicklly (00:06:21) in 407 generations to produce these parameter estimates :

Parameters = 
1.0e+02 * 
0.002865234375000 
8.700138633728027 
0.027473449707031 
0.009913467407227 
0.000003921508789 
0.002913574218750

with a fitness value :

fval = 
66.598426312801081

that is also the best of several runs, each with different 'InitialPopulationMatrix' amplitudes. The predicted ‘Y’ values are all complex.

Related Solutions

MATLAB: Can someone help me with these 3D plots

you can plot these shapes in the following two ways(you can choose better way you like):

%%Seashell

[u, v] = meshgrid(linspace(0, 6 * pi, 100), linspace(0, 2 * pi, 100));
x = 2 * (1 - exp(u / (6 * pi))) .* cos(u) .* (cos(0.5 * v))^2;
y = 2 * (-1 + exp(u / (6 * pi))) .* sin(u) .* (cos(0.5 * v))^2;
z = 1 - exp(u / (3 * pi)) - sin(v) + exp(u / (6 * pi)) .* sin(v);
figure, surf(x, y, z), xlabel('x'), ylabel('y'), zlabel('z'), title('Seashell')
%%Eight torus

c = 1;
[u, v] = meshgrid(linspace(-pi, pi, 100));
x = cos(u) .* (c + sin(v) .* cos(u) - sin(2 * v) .* (sin(u) / 2));
y = sin(u) .* (c + sin(v) .* cos(u) - sin(2 * v) .* (sin(u) / 2));
z = sin(u) .* sin(v) + cos(u) .* (sin(2 * v) / 2);
figure, surf(x, y, z), xlabel('x'), ylabel('y'), zlabel('z'), title('Eight torus')
%%Helical spring

r1 = 0.25;
r2 = 0.25;
T = 2;
n = 6;
[u, v] = meshgrid(linspace(0, 2 * n * pi, 100), linspace(0, 2 * pi, 100));
x = (1 - r1 * cos(v)) .* cos(u);
y = (1 - r1 * cos(v)) .* sin(u);
z = r2 * (sin(v) + T * u / pi);
figure, surf(x, y, z), xlabel('x'), ylabel('y'), zlabel('z'), title('Helical spring')
%%Cornucopia

a = 0.3;
b = 0.5;
[u, v] = meshgrid(linspace(0, 2 * pi, 100), linspace(-3, 3, 100));
x = exp(b * v) .* cos(v) + exp(a * v) .* cos(u) .* cos(v);
y = exp(b * v) .* sin(v) + exp(a * v) .* cos(u) .* sin(v);
z = exp(a * v) .* sin(u);
figure, surf(x, y, z), xlabel('x'), ylabel('y'), zlabel('z'), title('Cornucopia')

Or alternatively you can also plot them using ezsurf as follows:

syms u v
%%Seashell
x = 2 * (1 - exp(u / (6 * pi))) * cos(u) * (cos(0.5 * v))^2;
y = 2 * (-1 + exp(u / (6 * pi))) * sin(u) * (cos(0.5 * v))^2;
z = 1 - exp(u / (3 * pi)) - sin(v) + exp(u / (6 * pi)) * sin(v);
figure, ezsurf(x, y, z)
%%Eight torus
c = 1;
x = cos(u) * (c + sin(v) * cos(u) - sin(2 * v) * (sin(u) / 2));
y = sin(u) * (c + sin(v) * cos(u) - sin(2 * v) * (sin(u) / 2));
z = sin(u) * sin(v) + cos(u) * (sin(2 * v) / 2);
figure, ezsurf(x, y, z)
%%Helical spring
r1 = 0.25;
r2 = 0.25;
T = 2;
n = 6;
x = (1 - r1 * cos(v)) .* cos(u);
y = (1 - r1 * cos(v)) .* sin(u);
z = r2 * (sin(v) + T * u / pi);
figure, ezsurf(x, y, z)
%%Cornucopia
a = 0.3;
b =0.5;
x = exp(b * v) .* cos(v) + exp(a * v) .* cos(u) .* cos(v);
y = exp(b * v) .* sin(v) + exp(a * v) .* cos(u) .* sin(v);
z = exp(a * v) .* sin(u);
figure, ezsurf(x, y, z)

I hope it helps. Good luck!

MATLAB: Multiple variable analysis representation

[A, B, C, D] = ndgrid(1:10, 5:5:50, 0.1:0.1:1, 3:3:30);
y = A .* x.^3 + B .* x.^2 + C .* x + D;
[miny, minidx] = min(y(:));
[maxy, maxidx] = max(y(:));
parms_at_min = [A(minidx); B(minidx); C(minidx); D(minidx)];
parms_at_max = [A(maxidx); B(maxidx); C(maxidx); D(maxidx)];

If you have multiple X values then

numx = length(x);
[A, B, C, D, X] = ndgrid(1:10, 5:5:50, 0.1:0.1:1, 3:3:30, x);
y = A .* X.^3 + B .* X.^2 + C .* X + D;
rsy = reshape(y,[], numx);
[miny, minidx] = min(rsy);  %first dimension

[maxy, maxidx] = max(rsy);  %first dimension
parms_at_min = [A(minidx); B(minidx); C(minidx); D(minidx)];
parms_at_max = [A(maxidx), B(maxidx), C(maxidx), D(maxidx)];

I make no claims that this is the most efficient way.

Best Answer

Related Solutions

MATLAB: Can someone help me with these 3D plots

MATLAB: Multiple variable analysis representation

Related Question