Calculate mean and confidence interval
Show older comments
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.
Answers (1)
Walter Roberson
on 21 Jun 2021
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.
6 Comments
Rajvi Amle
on 21 Jun 2021
Edited: Rajvi Amle
on 21 Jun 2021
Walter Roberson
on 22 Jun 2021
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
on 24 Jun 2021
Rajvi Amle
on 27 Jun 2021
Walter Roberson
on 27 Jun 2021
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
on 27 Jun 2021
Edited: Rajvi Amle
on 27 Jun 2021
Categories
Find more on Loops and Conditional Statements in Help Center and File Exchange
on 21 Jun 2021
on 27 Jun 2021
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
