Why does lsqcurvefit result in complex parameters?

I want to fit experimental impedance data (real(impedance),imag(impedance),frequency).

The written function for getting the parameters has complex formulas for description of the impedance, but real-valued parameters. Unfortunately after the lsqcurvefit I get complex-valued parameters.

(The resnorm2 has a value of apprax 1e-6, which is not so bad.)

I want only real-valued parameters. What can I do to solve this problem?

My question is related to this question: https://ch.mathworks.com/matlabcentral/answers/320390-how-to-avoid-complex-eigenvalues-of-the-matrix-in-its-non-linear-regression-lsqcurvefit (but it was not helpful enough.) I've tried different things but without success till now and would be glad to get a solution from here.

---

In the code file "Fit3_exp.m" I wrote:

 % Start Parameters (initial guess)
 startPar = [0.45*1e-3, 1.4*1.0846e-7, 0.016*4.111, 0.7*1e-3, 11, 0.00020];
 % Make the best fitting
 fitPar = lsqcurvefit(@Fit3_f1, startPar, vertcat(f, f), vertcat(real(Z_exp), imag(Z_exp)),[0 0 0 0 0 0],[1e-3 1e6 1e6 1e6 1e6 1e6])
 % Vector of fitted Zfit data
 zFit3 = Fit3_f1(fitPar, vertcat(f, f));
 % Plot Fitted Line
 plot(zFit3(1:length(zFit3)/2),-zFit3(1+length(zFit3)/2:end), 'r','DisplayName','computed fitting data')

In the function file "Fit3_f1.m" I wrote:

function Z_fit3 = Fit3_f1(par,f)
file='G:\users\DA\Ueberblick EIS.xlsx';
exp_data_Freq = xlsread(file,'A4','A3:A73');
f = exp_data_Freq(:,1);
f = f(1:length(f)/1); % 
w = 2*pi*f ;

R0_fit0     = par(1);      
L_fit0      = par(2);    
R_bat_fit0  = par(3);      
R_ct_fit0   = par(4);   
C_dl_fit0   = par(5);      
s_w_fit0    = par(6);   

L_bat_fit0 = 1i*w*L_fit0;
R_el_fit0 = (L_bat_fit0*R_bat_fit0)./(L_bat_fit0+R_bat_fit0);
Z_d_fit0 = (s_w_fit0 ./ sqrt(w)) * (1-1i);        
Z_c_fit0 = -1i ./ (C_dl_fit0*w);                  
R_dl_fit0 =  ((Z_c_fit0 .* (R_ct_fit0 + Z_d_fit0)))./(R_ct_fit0 + Z_d_fit0 + Z_c_fit0);
Z_fit0 = R0_fit0 + R_el_fit0 + R_dl_fit0;             
Z_fit3 = vertcat(real(Z_fit0), imag(Z_fit0));

I get now only real-valued Parameters, however the fitting is really poor (resnorm2 = 3.5017e-04).

Are there mistakes in my formulas (e.g. with refering to the correct fit function (with indices) or is something else wrong)?

[Torsten]

Fit real and imaginary part of your data to real and imaginary part of your fit function. This way, you only work with real numbers, and the estimated parameters will remain real-valued.

[hamp]

I already tried it. Unfortunately, I don't get sense-making parameters. Please see the edited question.

[Torsten]

I don't see an obvious programming error, but you could test the fitting procedure by producing artificial data from your model, perturb these data slightly and try to reproduce the parameters you have used to generate these data.

If this works, you can be sure, that technically, your code is correct.

Then you can try "multistart" or modify your model.

Best wishes

Torsten.

[hamp]

It would be nice to find only real-valued parameter solutions. I'm still working on this issue.

[Torsten]

If you supply real-valued experimental data (real(Z_exp), imag(Z_exp)), return real-valued model data (real(Z_fit0), imag(Z_fit0)) and start with real-valued parameters, the parameters will remain real-valued throughout the fitting process.

Answer 1

[Walter Roberson]

It is because you are ignoring the second input parameter.

Note: for efficiency, do not read files in the objective function. Read the files beforehand and pass them into the objective function by creating an anonymous function that refers to them.

[hamp]

Thank you Walter for your input. I'm a beginner. Could you show me how the codes should look like?

[Walter Roberson]

I do not have your data to test with so I created some random data and changed to code to work with the complex values directly. The results I got back were real-valued for fitPar.

fit4_driver.m:

file='G:\users\DA\Ueberblick EIS.xlsx';
%f = xlsread(file,'A4','A3:A73');
f = rand(71,1);
%realpart = xlsread(file,'TS40','B3:B73');
%imagpart = xlsread(file,'TS40','C3:C73');
realpart = rand(71,1);
imagpart = randn(71,1);
Z_exp = complex(realpart, imagpart);

% Start Parameters (initial guess)
startPar = [0.45*1e-3, 1.4*1.0846e-7, 0.016*4.111, 0.7*1e-3, 11, 0.00020];
% Make the best fitting
fitPar = lsqcurvefit(@Fit4_f1, startPar, f, Z_exp, [0 0 0 0 0 0], [1e-3 1e6 1e6 1e6 1e6 1e6]);

% Vector of fitted Zfit data
zFit4 = Fit4_f1(fitPar, f);
% Plot Fitted Line
plot(real(zFit4),-imag(zFit4), 'r','DisplayName','computed fitting data')

Fit4_f1.m:

function Z_fit0 = Fit4_f1(par,f)
f = f(1:length(f)/1); % 
w = 2*pi*f ;

R0_fit0     = par(1);      
L_fit0      = par(2);    
R_bat_fit0  = par(3);      
R_ct_fit0   = par(4);   
C_dl_fit0   = par(5);      
s_w_fit0    = par(6);   

L_bat_fit0 = 1i*w*L_fit0;
R_el_fit0 = (L_bat_fit0*R_bat_fit0)./(L_bat_fit0+R_bat_fit0);
Z_d_fit0 = (s_w_fit0 ./ sqrt(w)) * (1-1i);        
Z_c_fit0 = -1i ./ (C_dl_fit0*w);                  
R_dl_fit0 =  ((Z_c_fit0 .* (R_ct_fit0 + Z_d_fit0)))./(R_ct_fit0 + Z_d_fit0 + Z_c_fit0);
Z_fit0 = R0_fit0 + R_el_fit0 + R_dl_fit0;             
%Z_fit4 = vertcat(real(Z_fit0), imag(Z_fit0));
end

Answer 2

[hamp]

I asked my colleague and was able to rewrite the code as shown below:

file='G:\users\DA\Ueberblick EIS.xlsx';
f = xlsread(file,'A4','A3:A73');
realpart = xlsread(file,'TS40','B3:B73');
imagpart = xlsread(file,'TS40','C3:C73');
Z_exp=realpart+i*imagpart;
realpart=real(Z_exp);
imagpart=imag(Z_exp);

% Start Parameters (initial guess)
startPar = [0.45*1e-3, 1.4*1.0846e-7, 0.016*4.111, 0.7*1e-3, 11, 0.00020];
% Make the best fitting
fitPar = lsqcurvefit(@Fit4_f1, startPar, vertcat(f, f), vertcat(real(Z_exp), imag(Z_exp)),[0 0 0 0 0 0],[1e-3 1e6 1e6 1e6 1e6 1e6])
% Vector of fitted Zfit data
zFit4 = Fit4_f1(fitPar, vertcat(f, f));
% Plot Fitted Line
plot(zFit4(1:length(zFit4)/2),-zFit4(1+length(zFit4)/2:end), 'r','DisplayName','computed fitting data')

function Z_fit4 = Fit4_f1(par,f)
f = f(1:length(f)/1); % 
w = 2*pi*f ;

R0_fit0     = par(1);      
L_fit0      = par(2);    
R_bat_fit0  = par(3);      
R_ct_fit0   = par(4);   
C_dl_fit0   = par(5);      
s_w_fit0    = par(6);   

L_bat_fit0 = 1i*w*L_fit0;
R_el_fit0 = (L_bat_fit0*R_bat_fit0)./(L_bat_fit0+R_bat_fit0);
Z_d_fit0 = (s_w_fit0 ./ sqrt(w)) * (1-1i);        
Z_c_fit0 = -1i ./ (C_dl_fit0*w);                  
R_dl_fit0 =  ((Z_c_fit0 .* (R_ct_fit0 + Z_d_fit0)))./(R_ct_fit0 + Z_d_fit0 + Z_c_fit0);
Z_fit0 = R0_fit0 + R_el_fit0 + R_dl_fit0;             
Z_fit4 = vertcat(real(Z_fit0), imag(Z_fit0));
end

I get following error code:

Error using lsqcurvefit (line 251)
Function value and YDATA sizes are not equal.

Error in Fit4_exp (line 12)
fitPar = lsqcurvefit(@Fit4_f1, startPar, vertcat(f, f), vertcat(real(Z_exp), imag(Z_exp)),[0 0 0 0 0 0],[1e-3 1e6 1e6 1e6 1e6 1e6])

length(f) = length(realpart) = length(imagpart) = 71. The dimensions of the Input values are the same.

What could I do to solve this issue?

[Walter Roberson]

The message is correct. You have

 Z_fit4 = vertcat(real(Z_fit0), imag(Z_fit0));

so you are producing a 284 x 1 output in response to a 142 x 1 input. Your fitting function must produce exactly one output for every input.

[hamp]

I believe, I don't understand how the function definitions: "function Z_fit4", "Fit4_f1", "Z_fit0" are used. Can you help me to understand?

[Walter Roberson]

function Z_fit4 = Fit4_f1(something)

means that the function name is Fit4_f1 and that whenever it is invoked, the value to be returned is whatever has been assigned to the variable Z_fit4. lsqcurvefit requires that the output be the same size as the input. Your current code is calculating an temporary variable Z_fit0 that is the right size, but splits it into real and imaginary parts and so is returning twice as much data as expected.

The current code is also not taking into account that half of the input corresponds to imaginary parts.

[hamp]

Thank you for your hints. I propose following code, which does result in real-valued parameters. What do you think about?:

clear
clc
close all

Z_real_exp = [3.1300e-03; 2.6849e-03; 2.2763e-03; 1.8945e-03; 1.6301e-03; 1.3839e-03; 1.1961e-03; 9.8979e-04; 8.4209e-04; 7.1154e-04; 6.2088e-04; 5.5614e-04; 5.0313e-04; 4.7850e-04; 4.6345e-04; 4.5015e-04; 4.5227e-04; 4.5836e-04; 4.6635e-04; 4.8297e-04; 4.9927e-04; 5.1492e-04; 5.5454e-04; 5.9414e-04; 6.4342e-04; 6.9217e-04; 7.4891e-04; 8.0603e-04; 8.6047e-04; 9.1073e-04; 9.5977e-04; 1.0020e-03; 1.0359e-03; 1.0676e-03; 1.0923e-03; 1.1135e-03; 1.1370e-03; 1.1582e-03; 1.1721e-03; 1.1896e-03; 1.2032e-03; 1.2162e-03; 1.2329e-03; 1.2466e-03; 1.2607e-03; 1.2782e-03; 1.2936e-03; 1.3068e-03; 1.3318e-03; 1.3512e-03; 1.3712e-03; 1.3991e-03; 1.4272e-03; 1.4509e-03; 1.4944e-03; 1.5179e-03; 1.5632e-03; 1.6116e-03; 1.6520e-03; 1.7130e-03; 1.7745e-03; 1.8599e-03; 1.9418e-03; 2.0433e-03; 2.1558e-03; 2.2881e-03; 2.4454e-03; 2.5902e-03; 2.7751e-03; 2.9336e-03; 3.1337e-03];
Z_imag_exp = [-9.3804e-03;   -7.7476e-03;   -6.5142e-03;   -5.3852e-03;   -4.5069e-03;   -3.7248e-03;   -3.0986e-03;   -2.5774e-03;   -2.1278e-03;   -1.7242e-03;   -1.4018e-03;   -1.1185e-03;   -8.7643e-04;   -6.9883e-04;   -5.3282e-04;   -3.8933e-04;   -2.8018e-04;   -1.8635e-04;   -1.0688e-04;   -3.7385e-05; 2.5667e-05; 7.7504e-05; 1.2548e-04; 1.7315e-04; 2.1478e-04; 2.4122e-04; 2.5791e-04; 2.6966e-04; 2.6735e-04; 2.5769e-04; 2.4846e-04; 2.3032e-04; 2.1241e-04; 1.9829e-04; 1.8456e-04; 1.6861e-04; 1.5832e-04; 1.5061e-04; 1.3993e-04; 1.3550e-04; 1.3094e-04; 1.2622e-04; 1.3010e-04; 1.3582e-04; 1.4180e-04; 1.4279e-04; 1.4525e-04; 1.5980e-04; 1.6961e-04; 1.8627e-04; 2.0286e-04; 2.2284e-04; 2.4247e-04; 2.7393e-04; 3.0023e-04; 3.3634e-04; 3.7935e-04; 4.2925e-04; 4.8946e-04; 5.5468e-04; 6.2873e-04; 7.2369e-04; 8.1695e-04; 9.3529e-04; 1.0378e-03; 1.1885e-03; 1.3078e-03; 1.4658e-03; 1.6081e-03; 1.7949e-03; 2.0005e-03];
Z_exp=Z_real_exp-1i*Z_imag_exp;

f=transpose(logspace(4,-3,71));

inputvariablevertcat = vertcat (f,f);
inputdatavertcat = vertcat(Z_real_exp,Z_imag_exp);

startPar = [0.45*1e-3, 1.4*1.0846e-7, 0.016*4.111, 0.7*1e-3, 11, 0.00020];
LB = [0, 0, 0, 0, 0, 0];
UB = [100, 100, 100, 100, 100, 100];
[fitPar,resnorm2,residual2,exitflag2,output2] = lsqcurvefit(@parestimation,startPar,inputvariablevertcat,inputdatavertcat,LB,UB) 

fitPar

function fitfunction = parestimation(par,inputvariablevertcat)

w = 2*pi*inputvariablevertcat ;
R0_fit0     = par(1);      
L_fit0      = par(2);    
R_bat_fit0  = par(3);      
R_ct_fit0   = par(4);   
C_dl_fit0   = par(5);      
s_w_fit0    = par(6);   
L_bat_fit0 = 1i*w*L_fit0;
R_el_fit0 = (L_bat_fit0*R_bat_fit0)./(L_bat_fit0+R_bat_fit0);
Z_d_fit0 = (s_w_fit0 ./ sqrt(w)) * (1-1i);        
Z_c_fit0 = -1i ./ (C_dl_fit0*w);                  
R_dl_fit0 =  ((Z_c_fit0 .* (R_ct_fit0 + Z_d_fit0)))./(R_ct_fit0 + Z_d_fit0 + Z_c_fit0);
Z_fit0 = R0_fit0 + R_el_fit0 + R_dl_fit0;

firsthalf = Z_fit0(1:(length(Z_fit0)/2));
secondhalf = Z_fit0((length(Z_fit0)/2+1):length(Z_fit0));

fitfunction = vertcat(real(firsthalf),imag(secondhalf)); 

end