Calculate mean and confidence interval

I am trying to calculate mean and confidence interval of x=1x10 cells containing 146x18 simulations calculated in a loop for 10 iterations. Although when code is run, it shows me the error that

'matrix dimensions must agree'.

Error in toy_example (line 262)

CI95 = mean_x + ts.*(std_x/sqrt(length(x)));

Could someone please help me with this suggesting what mistake I did. This would be really very helpful and appreciated.

My code for the loop for 10 iterations and them mean and CI for those 10 simulations/iterations is mentioned below and also attached the values of 'x':

%Enteries in a loop
n=10;
t=[];
x=[];

for i=1:n
    tspan=[0 max(EXP.t)]; 
    MAT=[SYSTEM.m_BW(i) SYSTEM.w_L(i) SYSTEM.w_K(i) SYSTEM.w_Lu(i) SYSTEM.w_BR(i) SYSTEM.w_Blood(i) SYSTEM.w_Plasma(i) SYSTEM.m_L(i) SYSTEM.m_BR(i) SYSTEM.m_K(i) SYSTEM.m_S(i) SYSTEM.m_Lu(i) SYSTEM.m_Blood(i) SYSTEM.m_Plasma(i)...
        SYSTEM.V_L(i) SYSTEM.V_BR(i) SYSTEM.V_K(i) SYSTEM.V_S(i) SYSTEM.V_Lu(i) SYSTEM.V_Blood(i) SYSTEM.V_Plasma(i) SYSTEM.F_C(i) SYSTEM.F_L(i) SYSTEM.F_BR(i) SYSTEM.F_K(i) SYSTEM.F_S(i) SYSTEM.F_Res(i) SYSTEM.F_bal(i) SYSTEM.F_Bile(i) SYSTEM.F_Urine(i)...
        SYSTEM.V_Res(i) SYSTEM.V_bal(i) SYSTEM.V_L_b(i) SYSTEM.V_L_t(i) SYSTEM.V_BR_b(i) SYSTEM.V_BR_t(i) SYSTEM.V_K_b(i)  SYSTEM.V_K_t(i) SYSTEM.V_S_b(i) SYSTEM.V_S_t(i) SYSTEM.V_Lu_b(i) SYSTEM.V_Lu_t(i) SYSTEM.V_Res_b(i) SYSTEM.V_Res_t(i)... 
        DRUG.m_Au_iv(i) DRUG.M_Au_iv(i)]';
     options=odeset('AbsTol',10e-2,'RelTol',10e-2,'Stats','on');
     col=zeros(18,1);
     tic;
     [t0,x0]=ode15s(@ode_toy, tspan, col,[],ov,DRUG,MAT);    
       t=[t, {t0}];
       x=[x, {x0}]; 
      toc; 
        
end

% Confidence interval of 95% for simulated data

mean_x=mean([x{:,2}]);               % convert cell array to matrix, mean by row
% N=size(y,1);
std_x=mean_x/size(x,2);                     % std dev of mean(y) is mean(y)/nObs
SEM=(std_x/sqrt(length(x)))
ts = tinv([0.025  0.975],length(x)-1);      % T-Score
CI95 = mean_x + ts.*(std_x/sqrt(length(x)));

Please do help me.

Answer 1

[Walter Roberson]

mean_x=mean([x{:,2}]);               % convert cell array to matrix, mean by row

Your x is a cell array, with entries that might be different number of rows because the ode*() routines can output a variable number of rows when you provide a time span that is exactly two elements.

The syntax you are using [x{:,2}] is going to access the second element of the cell array only, and extra the one matrix there. Your initial conditions are length 18, so that is going to be a something-by-18 matrix, Then the [] around the single something-by-18 matrix is going to leave it as a something-by-18 matrix. Your x{:,2} looks on the surface like it might be selecting multiple entries in x, but x is a row vector of cell arrays because of the way you build it with [x,{x0}] so there is only one row in x and so x{:,2} is going to be x{1,2}

You then take the mean of that something-by-18 matrix along the first non-scalar dimension, so you will get a 1 x 18 vector.

std_x=mean_x/size(x,2);                     % std dev of mean(y) is mean(y)/nObs

size(x,2) is going to be the number of cells in the cell array, which will be n. But why would you divide the mean created from a single cell entry by the number of iterations you did?

SEM=(std_x/sqrt(length(x)))

why are you mixxing using length(x) and size(x,2) to get the number of entries in x ? Be consistent or people will get stuck trying to figure out if there are any circumstances under which the two could differ.

ts = tinv([0.025  0.975],length(x)-1);      % T-Score

That is going to be a vector of length 2.

CI95 = mean_x + ts.*(std_x/sqrt(length(x)));

mean_x is 1 x 18, std_x is 1 x 18. sqrt(length(x)) is scalar, so you have 1 x 18 + 1 x 2 .* 1 x 18 . But 1 x 2 .* 1 x 18 is going to fail.

If you made

ts = tinv([0.025; 0.975],length(x)-1);      % T-Score

then it would be 2 x 1 and that can be .* with 1 x 18 giving 2 x 18, which can be added to 1 x 18, giving 2 x 18.

[Rajvi Amle]

Hello @Walter Roberson. mean_x=mean([x{:,2}]); isn't working and plot of CI shows straight diagnol line which isnt the actual result. So I made some changes in the tried to calculate mean using cell2mat.

%Enteries in a loop
n=10;
t=[];
x=[];
for i=1:n
    tspan=[0 max(EXP.t)]; 
    MAT=[SYSTEM.m_BW(i) SYSTEM.w_L(i) SYSTEM.w_K(i) SYSTEM.w_Lu(i) SYSTEM.w_BR(i) SYSTEM.w_Blood(i) SYSTEM.w_Plasma(i) SYSTEM.m_L(i) SYSTEM.m_BR(i) SYSTEM.m_K(i) SYSTEM.m_S(i) SYSTEM.m_Lu(i) SYSTEM.m_Blood(i) SYSTEM.m_Plasma(i)...
        SYSTEM.V_L(i) SYSTEM.V_BR(i) SYSTEM.V_K(i) SYSTEM.V_S(i) SYSTEM.V_Lu(i) SYSTEM.V_Blood(i) SYSTEM.V_Plasma(i) SYSTEM.F_C(i) SYSTEM.F_L(i) SYSTEM.F_BR(i) SYSTEM.F_K(i) SYSTEM.F_S(i) SYSTEM.F_Res(i) SYSTEM.F_bal(i) SYSTEM.F_Bile(i) SYSTEM.F_Urine(i)...
        SYSTEM.V_Res(i) SYSTEM.V_bal(i) SYSTEM.V_L_b(i) SYSTEM.V_L_t(i) SYSTEM.V_BR_b(i) SYSTEM.V_BR_t(i) SYSTEM.V_K_b(i)  SYSTEM.V_K_t(i) SYSTEM.V_S_b(i) SYSTEM.V_S_t(i) SYSTEM.V_Lu_b(i) SYSTEM.V_Lu_t(i) SYSTEM.V_Res_b(i) SYSTEM.V_Res_t(i)... 
        DRUG.m_Au_iv(i) DRUG.M_Au_iv(i)]';
     options=odeset('AbsTol',10e-2,'RelTol',10e-2,'Stats','on');
     col=zeros(18,1);
     tic;
     [t0,x0]=ode15s(@ode_toy, tspan, col,[],ov,DRUG,MAT);    
       t=[t, {t0}];
       x=[x, {x0}];
      toc; 
        
end

% Confidence interval of 95% for simulated data

mean_x=mean(cell2mat(x),2);               % convert cell array to matrix, mean by row
% N=size(y,1);
std_x=mean_x/size(x,2);                     % std dev of mean(y) is mean(y)/nObs
ts = tinv([0.025;  0.975],length(x)-1);      % T-Score
CI95 = mean_x + ts.*(std_x/sqrt(length(x)));


% 
figure
plot(mean_x,'-b')                                  % Plot Mean Of All Experiments
hold on
plot(CI95,'--b')                                    % Plot 95% Confidence Intervals Of All Experiments
hold off
grid
legend('\mu','95% Confidence Limits','Location','best')

when run, shows the error as:

Error using cat
Dimensions of arrays being concatenated are not consistent.

Error in cell2mat (line 75)
            m{n} = cat(2,c{n,:});

Error in toy_example (line 257)
mean_x=mean(cell2mat(x),2);               % convert cell array to matrix, mean by row

What can be the problem?

[Walter Roberson]

As I wrote:

Your x is a cell array, with entries that might be different number of rows because the ode*() routines can output a variable number of rows when you provide a time span that is exactly two elements.

The implicit hint there is that you should instead provide a vector for tspan that has at least three elements. ode*() will output one row for each of the times given in the vector, and so it would be a fixed number of rows.

I am not at all clear what you are trying to take the mean() over. If you put all of the x together into an array with 12 columns, then mean() along the first dimension would give you a 1 x 12 result that was the mean over all of the times over all of the runs, which does not sound to me to be beneficial. If you put all of the x together into an array that has (n * 12) columns, then mean() along the first dimension would give you a 1 x (n*12) result that was the mean over all of the times, but broken up by runs... taking the mean over all the time entries for one run does not sound likely to be beneficial to me.

If you put all of the x together into an array with 12 columns, then mean() along the second dimension, you would be taking a mean along the different variables of integration, which does not sound likely to be beneficial to me.

If you put all of the x together into an array with (n*12) columns, then mean() along the second dimension, you would be taking a mean along the different variables of integration over all of the n, which does not sound likely to be beneficial to me.

I could see possibilities if you were taking the final row of each of the n and taking the mean over the n, which would have possibilities to explore the confidence intervals about where each of the boundary conditions ends up, subject to uncertainties introduced through MAT values.

[Rajvi Amle]

@Walter Roberson I understood what you suggested. I am calculating mean of values of x(:,1)+x(:,2) for every iterations from x as x contains 25x18 simulation results over t with 25x1 matrix. So I want values x(:,1) to be converted to values in row and then calculate its mean. Similarly for x(:,2). I tried in two ways to get my mean.

1st way

t=cell(1,10); y=t;                                        % the output cell arrays as row 
figure
hold on

for i=1:n
    tspan=[0:24]; 
    MAT=[SYSTEM.m_BW(i) SYSTEM.w_L(i) SYSTEM.w_K(i) SYSTEM.w_Lu(i) SYSTEM.w_BR(i) SYSTEM.w_Blood(i) SYSTEM.w_Plasma(i) SYSTEM.m_L(i) SYSTEM.m_BR(i) SYSTEM.m_K(i) SYSTEM.m_S(i) SYSTEM.m_Lu(i) SYSTEM.m_Blood(i) SYSTEM.m_Plasma(i)...
        SYSTEM.V_L(i) SYSTEM.V_BR(i) SYSTEM.V_K(i) SYSTEM.V_S(i) SYSTEM.V_Lu(i) SYSTEM.V_Blood(i) SYSTEM.V_Plasma(i) SYSTEM.F_C(i) SYSTEM.F_L(i) SYSTEM.F_BR(i) SYSTEM.F_K(i) SYSTEM.F_S(i) SYSTEM.F_Res(i) SYSTEM.F_bal(i) SYSTEM.F_Bile(i) SYSTEM.F_Urine(i)...
        SYSTEM.V_Res(i) SYSTEM.V_bal(i) SYSTEM.V_L_b(i) SYSTEM.V_L_t(i) SYSTEM.V_BR_b(i) SYSTEM.V_BR_t(i) SYSTEM.V_K_b(i)  SYSTEM.V_K_t(i) SYSTEM.V_S_b(i) SYSTEM.V_S_t(i) SYSTEM.V_Lu_b(i) SYSTEM.V_Lu_t(i) SYSTEM.V_Res_b(i) SYSTEM.V_Res_t(i)... 
        DRUG.m_Au_iv(i) DRUG.M_Au_iv(i)]';
     options=odeset('AbsTol',10e-2,'RelTol',10e-2,'Stats','on');
%      col=zeros(18,1);
     x0=zeros(18,1);
     tic;
%      [t0,x0]=ode15s(@ode_toy, tspan, col,[],ov,DRUG,MAT);    
    [t{i},x]=ode15s(@ode_toy, tspan, x0,[],ov,DRUG,MAT); 
%        t=[t, {t0}];
%        x=[x, {x0}]; 
       y{i}=(x(:,2)+x(:,1));
    plot(t{i},y{i},'-.',EXP.t,EXP.Blood_15,'o','Color', 'blue', 'LineWidth', 2)
    lgd= legend('simulated data','experimental data','Location','northeast');
    ylabel(['Amount of GSH-AuNPs in blood [mg/g]'])
    xlabel(['Time [h]'])
      toc; 
        
end


% Confidence interval of 95% for simulated data

mean_y=mean(cell2mat(y),2);               % convert cell array to matrix, mean by row
% N=size(y,1);
std_y=mean_y/size(y,2);                     % std dev of mean(y) is mean(y)/nObs
ts = tinv([0.025  0.975],length(y)-1);      % T-Score
CI95 = mean_y + ts.*(std_y/sqrt(length(y)));
% 
% 
% % 
figure
plot(mean_y,'-b')                                   % Plot Mean Of All Experiments
hold on
plot(CI95,'--b')                                    % Plot 95% Confidence Intervals Of All Experiments
hold off
grid
legend('mean','95% Confidence Limits','Location','best')
ylabel(['Amount of GSH-AuNPs in blood [mg/g]'])
xlabel(['Time [h]'])

which plots not the favourable result for what I want and 2nd way was:

n=10;
t=[];
x=[];

for i=1:n
    tspan=[0 max(EXP.t)]; 
    MAT=[SYSTEM.m_BW(i) SYSTEM.w_L(i) SYSTEM.w_K(i) SYSTEM.w_Lu(i) SYSTEM.w_BR(i) SYSTEM.w_Blood(i) SYSTEM.w_Plasma(i) SYSTEM.m_L(i) SYSTEM.m_BR(i) SYSTEM.m_K(i) SYSTEM.m_S(i) SYSTEM.m_Lu(i) SYSTEM.m_Blood(i) SYSTEM.m_Plasma(i)...
        SYSTEM.V_L(i) SYSTEM.V_BR(i) SYSTEM.V_K(i) SYSTEM.V_S(i) SYSTEM.V_Lu(i) SYSTEM.V_Blood(i) SYSTEM.V_Plasma(i) SYSTEM.F_C(i) SYSTEM.F_L(i) SYSTEM.F_BR(i) SYSTEM.F_K(i) SYSTEM.F_S(i) SYSTEM.F_Res(i) SYSTEM.F_bal(i) SYSTEM.F_Bile(i) SYSTEM.F_Urine(i)...
        SYSTEM.V_Res(i) SYSTEM.V_bal(i) SYSTEM.V_L_b(i) SYSTEM.V_L_t(i) SYSTEM.V_BR_b(i) SYSTEM.V_BR_t(i) SYSTEM.V_K_b(i)  SYSTEM.V_K_t(i) SYSTEM.V_S_b(i) SYSTEM.V_S_t(i) SYSTEM.V_Lu_b(i) SYSTEM.V_Lu_t(i) SYSTEM.V_Res_b(i) SYSTEM.V_Res_t(i)... 
        DRUG.m_Au_iv(i) DRUG.M_Au_iv(i)]';
     options=odeset('AbsTol',10e-2,'RelTol',10e-2,'Stats','on');
     col=zeros(18,1);
     tic;
     [t0,x0]=ode15s(@ode_toy, tspan, col,[],ov,DRUG,MAT);     
       t=[t, {t0}];
       x=[x, {x0}]; 
      toc; 
        
end


% Confidence interval of 95% for simulated data

mean_y=mean([x{:,2}]);                % convert cell array to matrix, mean by row
% N=size(y,1);
std_y=mean_y/size(x,2);                     % std dev of mean(y) is mean(y)/nObs
ts = tinv([0.025;  0.975],length(x)-1);      % T-Score
CI95 = mean_y + ts.*(std_y/sqrt(length(x)));
% 
% 
% % 
figure
plot(mean_y,'-b')                                   % Plot Mean Of All Experiments
hold on
plot(CI95,'--b')                                    % Plot 95% Confidence Intervals Of All Experiments
hold off
grid
legend('mean','95% Confidence Limits','Location','best')
ylabel(['Amount of GSH-AuNPs in blood [mg/g]'])
xlabel(['Time [h]'])

Also does not show any results.

The loop shouls show the graph for 10 iterations as illustrated in (Fig6.png) along with the experimental values and mean and confidence interval.

To make you understand better, I have attached an image of how should the mean be calculated in written.

Please do help me and suggest what I can do to calculate mean.

[Rajvi Amle]

@Walter Roberson Looking forward to your positive response. Could you please help me with the above question?

[Walter Roberson]

Please post your ode function. It looks unlikely at the moment that what you are asking for is physically meaningful, but not completely impossible.

[Rajvi Amle]

@Walter Roberson I have attached 3 files. 1st file toy_example have all defined variables and 2 biological variables calculated with random function, a loop for 10 iterations, post processing sections that are stored as results which can be used to plot as simulation results and in the end are given all the plottings. 2nd file ode_toy has all 18 ode equations and 3rd file is ae_toy which has algebraic equations and connected with the ODEs. All these three files are interconnected to each other. My model is quite complex but please do help me if possible. I only have problem with the calculation of mean as I am facing problem with the matrix size, dimensions, cell arrays etc.