How to graph a convolution with same time length as inputs

Question

Lilly 2023-10-11

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2032159-how-to-graph-a-convolution-with-same-time-length-as-inputs

评论： Lilly 2023-10-12

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() ?

Lilly 2023-10-11

Hello, yes I have but the graph when using 'same' doesn't look the way that it should.

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Matt J 2023-10-11

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2032159-how-to-graph-a-convolution-with-same-time-length-as-inputs#answer_1331199

编辑：Matt J 2023-10-11

在 MATLAB Online 中打开

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 个评论
显示 3更早的评论隐藏 3更早的评论

Paul 2023-10-12

在 MATLAB Online 中打开

But when using 'same', the values in Ci don't correspond to the values in tnew. In the plot below, Cisame has to be shifted to the right by 295 seconds to be consisten with time reference for the signals being convolved.

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

Cisame = conv(cpnew,vnew,'same')*dt;

Cifull = conv(cpnew,vnew,'full')*dt;

plot(tnew,cpnew,tnew,Cisame,(0:numel(Cifull)-1)*dt,Cifull);

title("Concentration over Time");

xlabel('Time (mins)');

ylabel('Concentration')

legend("Cp(t)","Cisame(t)","Cifull(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

Matt J 2023-10-12

编辑：Matt J 2023-10-12

Well, v is not a signal, but rather, a convolution kernel so it doesn't really "start" at a particular time. If we assume it is the impulse response of a causal system then, yes, we would normally assume everything starts at t=0.

Lilly 2023-10-12