- 'median' function: https://in.mathworks.com/help/matlab/ref/median.html
- 'bootci' function: https://in.mathworks.com/help/stats/bootci.html
How to plot median Confidence interval and how to comape two model ? im trynig to find the diffrence bw 2 models under difrent error. Here error are diff distbution,
12 Ansichten (letzte 30 Tage)
Ältere Kommentare anzeigen
N= 10000; % Number of observation
jN = 20; % Number of JumpsF
jumpInd = randperm(N,jN);
jumpInd = sort(jumpInd,'ascend');
jumpHeights = random('logn',log(15),0.5,jN+1,1);
jumpInd = [0 jumpInd N];
f = zeros(N,1);
a =[0.1 0.2 0.1 1.05]; % Real Values
k = a(4);
ind=zeros(N,1);
interval_N = N/5; %requires N to be divisible by 5
S=10; % No of simulation
alpha = 0.05;
for i=1:5
ind(interval_N*(i-1)+1:interval_N*(i-1)+interval_N/2)=1;
ind=logical(ind);
c = zeros(S,4);
end
for j= 1:S
for i = 1:jN+1
f(jumpInd(i)+1:jumpInd(i+1)) = jumpHeights(i);
end
v1=a(1)+a(2)*f+a(3)*f.^2; % Model One
v2=k*(v1); % Model Two
v3 = zeros(size(v1));
v3(ind)= a(1) + a(2)*f(ind) + a(3)*f(ind).^2;
v3(~ind)=k*(a(1) + a(2)*f(~ind) + a(3)*f(~ind).^2);
%v3= v3 + 7*randn(size(v3));
v3= v3 + 3*random('t',4,size(v3)); % T distribution
%v3 = v3 + 5* random('ev',0,1,size(v3)); % extream values distribution
%v3 = v3 + 2*poissrnd(v3);
% Defines the residual function.
% Here we let a1 = c(1), a2 = c(2), a3 = c(3) and k = c(4)
F = @(c) ind.*(c(1) + c(2)*f + c(3)*f.^2 - v3) +(1-ind).*(c(4)*(c(1) + ...
c(2)*f + c(3)*f.^2) - v3);
% Here we set options for lsqnonlin
options = optimoptions(@lsqnonlin,'Algorithm','trust-region-reflective');
% Here we set the initial guess for the parameters
c0 = [0.1,0.2,0.1 1]; % Initail values
[c(j,:),resnorm,residual,output] = lsqnonlin(F,c0,[-1e6,-1e6,-1e6,-1e6],[1e10,1e10,1e10,2],options);
% To check the uncertanity in the estimated parameters
v3_fit = zeros(size(v1));
v3_fit(ind)= c(j,1) + c(j,2)*f(ind) + c(j,3)*f(ind).^2;
v3_fit(~ind)=c(j,4)*(c(j,1) + c(j,2)*f(~ind) + c(j,3)*f(~ind).^2);
end
mu= mean(c(:,4));
sigma = std(c(:,4));
Median = median(c(:,4));
ci95 = cint_median(c(:,4),alpha);
e = v3 -v3_fit;
table(mu,sigma,Median,ci95)
%%%%%%%%%%%%%%%%%%% function to determine the confidence interval of median%%%%%%%%%%%%%%%%
function [interval] = cint_median(M,alpha)
%creates a confidence interval for the median of a dataset using a sign
%test based on binomial distribution.
if nargin < 2
alpha = 0.05;
end
%sign test
n = length(M);
m = median(M);
%confidence intervall based on median
c = icdf('bino', alpha/2 , n , 0.5);
s = sort(M);
interval = [s(c) m s(end - c)];
end
%%%%%%%%% least square method for model1%%%%%
function coefficient = LSMcurvefit_v1(x,y)
A = [ones(size(x)) x x.^2];
coefficient=(A'*A)\(A'*y);
end,
%%%%%%%%%%%%%%%Least square method for model 2%%%%%
function coefficient = LSMcurvefit_v2(x,y,c)
z = c(1)+c(2)*x+c(3)*x.^2;
coefficient = sum(z.*y)/sum(z.^2);
end
0 Kommentare
Antworten (1)
Balavignesh
am 9 Okt. 2023
Hi Muhammad,
As per my understanding, you would like to plot median confidence interval of the models.
I would suggest you to calculate the 'median' of your data using 'median' function. You could calculate the confidence interval using the 'bootci' or 'prctile' function to calculate the confidence interval. You could use 'patch' or 'fill' function to plot the confidence interval as a shaded area around the median line.
Here is an example code that could help you understand this:
% This is an example datast. Use your own dataset
load carsmall;
weights = Weight
% Compute the median and confidence interval
medianValue = median(weights);
ci = bootci(1000, @median, weights);
% Plotting the Median Confidence Interval
figure;
plot(1, medianValue, 'ko', 'MarkerSize', 8, 'LineWidth', 2);
hold on;
patch([0.9, 1.1, 1.1, 0.9], [ci(1), ci(1), ci(2), ci(2)], 'b', 'FaceAlpha', 0.3);
% Add labels and customize the plot
xlabel('Car');
ylabel('Weight (lbs)');
title('Median Weight with Confidence Interval');
% Remove the extra tick mark
set(gca, 'XTick', 1);
% Remove the box around the plot
box off;
hold off;
Please refer to the following documentation links to have more information on:
Hope that helps!
Balavignesh.
0 Kommentare
Siehe auch
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!