计算不同参数值的模拟平均值 [英] compute the average over simulation for different parameters values

问题描述

我正在研究类似于以下示例的内容: 我想计算<x(t)>，这是模拟次数中函数x(t)的平均值.为此，我生成以下代码:

I am working on something similar to the following example: I want to compute <x(t)> that is the average of a function x(t) over number of simulations. To do this, I generate the following code:

sim=50;% number of simulations 
t=linspace(0,1);% time interval

a_range=[1,2,3];% different values for the parameter a
b_range=[0,0.5,1];% different values for the parameter b
z=zeros(1,sim);
theta=zeros(1,sim);

for nplot=1:3
   a=a_range(nplot);
   b=b_range(nplot);
   average_x=zeros(nplot,sim);
   for i=1:sim
       z(i)=rand(1);% random number for every simulation
       theta(i)=pi*rand(1);% random number for every simulation    
    x=z(i)*t.^2+a*sin(theta(i))+b.*tan(theta(i));% the function

   end
   average_x(nplot,sim)=mean(x);% average over the number of simulations
end

fname=['xsin.mat'];
save(fname)

时间是1乘以100的向量，x是1乘以100的向量，average_x是1乘50的向量.我要寻找的是编写一个脚本来加载文件并针对时间对不同参数绘制平均值a和b.所以我想编写一个代码来生成三个数字，以便在图1中绘制平均值

The time is a vector 1 by 100 and x is a vector 1 by 100, and average_x is 1 by 50. What I am looking for is to write a script to load the file and plot the average against time for different parameters a and b. So I want to write a code to generate three figures such that in figure 1 I will plot the average

plot(t,average_x)

对于a = 1和b = 0.

for a=1 and b=0.

然后在图2中，我将再次绘制平均值，但对于a = 2和b = 0.5，依此类推.问题在于时间t的维数与平均值不相同.如何解决此问题并生成三个不同的数字.

Then in figure 2 I will plot the average again but for a=2 and b=0.5 and so on. The problem is the dimension of time t and the average are not the same. How can I fix this problem and generate three distinct figures.

推荐答案

如果我正确理解了您的意图，这就是您要寻找的:

If I got your intention correctly, this is what you look for:

sim = 50;% number of simulations 
t = linspace(0,1);% time interval

a_range = [1,2,3];% different values for the parameter a
b_range = [0,0.5,1];% different values for the parameter b

% NO NEED TO GENERATE THE RANDOM NUMBERS ONE BY ONE:
theta = pi*rand(sim,1);% random number for every simulation
z = rand(sim,1); % random number for every simulation

% YOU SOULD INITIALIZE ALL YOUR VARIABLES OUTSIDE THE LOOPS:
x = zeros(sim,numel(t)); 
average_x = zeros(3,numel(t));% the mean accross simulations
% for average accros time use:
% average_x = zeros(3,sim);

for nplot=1:3
   a = a_range(nplot);
   b = b_range(nplot);
   for i=1:sim
       x(i,:) = z(i)*t.^2+a*sin(theta(i))+b.*tan(theta(i));% the function
   end
   average_x(nplot,:) = mean(x); % average over the number of simulations
   % average_x(nplot,:) = mean(x,2); % average accross time
end
% save the relevant variables:
save('results.mat','average_x','t')

在另一个文件中，您可以编写:

In another file you can write:

load('results.mat')
for k = 1:size(average_x,1)
    figure(k)
    plot(t,average_x(k,:))
    title(['Parameter set ' num2str(k)])
    xlabel('Time')
    ylabel('mean x')
end

这是一个图中的图(如果您想在模拟中取平均值):

This is the plot in one figure (if you want then average over simulations):

顺便说一句，如果您想使代码更紧凑，更快速，则可以对其进行矢量化，主要是使用 bsxfun .这是您的代码的演示:

BTW, if you want to make your code more compact and fast, you can vectorize it, mainly using bsxfun. Here is a demonstration with your code:

% assuming all parameters are defined as above:
zt = bsxfun(@times,z,t.^2); % first part of the function 'z(i)*t.^2'
% second part of the function 'a*sin(theta(i)) + b.*tan(theta(i))':
ab = bsxfun(@times,a_range,sin(theta)) + bsxfun(@times,b_range,tan(theta));
% convert the second part to the right dimensions and size:
ab = repmat(reshape(ab,[],1,3),1,numel(t),1); 
x = bsxfun(@plus,zt,ab); % the function
average_x = squeeze(mean(x)); % take the mean by simulation
plot(t,average_x) % plot it all at once, as in the figure above
xlabel('Time')
ylabel('mean x')

这篇关于计算不同参数值的模拟平均值的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！