How to graph a convolution with same time length as inputs
2 次查看(过去 30 天)
显示 更早的评论
Hello,
So far this is my code. What I'm trying to do is plot Ci and Cp over time, but Ci is a convolution of Cp and v. Both Cp and v are of length 27, so Ci is of length 53 since it is a convolution. The issue I'm having here is that I don't know how to get an accurate representation of Ci when graphing since time is given by specific points (see t). These specific time points are also directly related to Cp. What's the best way to do this?
Thank you in advance!
%v denotes the exponential function
v=myfunction(t);
%conce is the function for Cp(t)
conce=Cp(t);
%CC is convolution of Cp(t)*the exponential function
CC=conv(conce,v);
%Ci is the convolution having the same size of the two input vectors
Ci=CC(1:27);
plot(t,Ci,t,conce);
title("Concentration over Time");
xlabel('Time (mins)');
ylabel('Concentration')
legend("Ci(t)","Cp(t)");
%Based on plot, it appears that the highest signal strength in the brain is
%around 39.93 minutes, which should be the best time to run the PET scan.
function conce=Cp(a)
t=[0,1.08,1.78,2.30,2.75,3.30,3.82,4.32,4.80,5.28,5.95,6.32,6.98,9.83,16.30,20.25,29.67,39.93,58,74,94,100,200,300,400,500,591];
conc=[0,84.9,230,233,220,236.4,245.1,230.0,227.8,261.9,311.7,321,316.6,220.7,231.7,199.4,211.1,190.8,155.2,140.1,144.2,139.737,111.3006,82.8375,54.3744,25.9113,0.0098];
conce=conc(t==a);
end
function v=myfunction(t);
k1=0.102;
k2=0.130;
k3=0.062;
k4=0.0068;
alpha1=(k2+k3+k4+(sqrt((k2+k3+k4)^2-(4*k2*k4))))/2;
alpha2=(k2+k3+k4-(sqrt((k2+k3+k4)^2-(4*k2*k4))))/2;
A=(k1*(k4+k3-alpha1))/(alpha2-alpha1);
B=(k1*(alpha2-k4-k3))/(alpha2-alpha1);
v=A*exp(-alpha1*t)+B*exp(-alpha2*t);
%Ci=Cp(t).*(A.*exp(-1*alpha1*t)+B.*exp(-1*alpha2*t));
end
2 个评论
Walter Roberson
2023-10-11
Have you considered looking at the 'same' and 'valid' options of conv() ?
采纳的回答
Matt J
2023-10-11
编辑:Matt J
2023-10-11
Convolution doesn't make sense if your time points are not equi-spaced. You should interpolate everything so that they are:
t=[0,1.08,1.78,2.30,2.75,3.30,3.82,4.32,4.80,5.28,5.95,6.32,6.98,9.83,16.30,20.25,29.67,39.93,58,74,94,100,200,300,400,500,591];
%v denotes the exponential function
v=myfunction(t);
%conce is the function for Cp(t)
cp=Cp(t);
dt=min(diff(t));
tnew=min(t):dt:max(t);
vnew=interp1(t,v,tnew);
cpnew=interp1(t,cp,tnew);
%CC is convolution of Cp(t)*the exponential function
Ci=conv(cpnew,vnew,'same')*dt;
plot(tnew,cpnew,tnew,Ci);
title("Concentration over Time");
xlabel('Time (mins)');
ylabel('Concentration')
legend("Cp(t)","Ci(t)");
%Based on plot, it appears that the highest signal strength in the brain is
%around 39.93 minutes, which should be the best time to run the PET scan.
function conce=Cp(a)
t=[0,1.08,1.78,2.30,2.75,3.30,3.82,4.32,4.80,5.28,5.95,6.32,6.98,9.83,16.30,20.25,29.67,39.93,58,74,94,100,200,300,400,500,591];
conc=[0,84.9,230,233,220,236.4,245.1,230.0,227.8,261.9,311.7,321,316.6,220.7,231.7,199.4,211.1,190.8,155.2,140.1,144.2,139.737,111.3006,82.8375,54.3744,25.9113,0.0098];
conce=conc(t==a);
end
function v=myfunction(t);
k1=0.102;
k2=0.130;
k3=0.062;
k4=0.0068;
alpha1=(k2+k3+k4+(sqrt((k2+k3+k4)^2-(4*k2*k4))))/2;
alpha2=(k2+k3+k4-(sqrt((k2+k3+k4)^2-(4*k2*k4))))/2;
A=(k1*(k4+k3-alpha1))/(alpha2-alpha1);
B=(k1*(alpha2-k4-k3))/(alpha2-alpha1);
v=A*exp(-alpha1*t)+B*exp(-alpha2*t);
%Ci=Cp(t).*(A.*exp(-1*alpha1*t)+B.*exp(-1*alpha2*t));
end
5 个评论
更多回答(0 个)
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Line Plots 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!