Help with basic code

1 次查看（过去 30 天）
Susan
Susan 2011-3-22
I attached the pdf its the last page 9 I meant to code that in MATLAB, I believe even the code is basic but I am not sure why is not working. I tried to change the code several time but no luck.
http://www.stat.sdu.dk/matstat/yuri-st505/w5slides1.pdf
from the Paper:
(1) t = 0, I = 0.
(2) Generate a random number U ~ U(0, 1).
(3) t = t - (1/Lam)*ln(U). If t > T then stop.
(4) Generate a random number U ~ U(0, 1).
(5) If U <= Lam(t)/Lam, set I = I + 1, S(I) = t.
(6) Go to step 2.
Output:
I: the number of events at time T ,
S(1), . . . , S(I): the event times.
2. Direct generation of successive event times
function nonhomogeneous (i,s)
i = 0;
t = 010;
T = length(cos(t));
s = [0,t];
u = rand (1,1);
t = t - (1 /cos(t)*log(u));
while t <= T
u2 = rand (1,1);
if ( u2 <= cos(t) / max(cos(t)))
hold on
plot(t,cos(t))
i = i+1;
s(i)=t;
u = rand (1,1);
t = t - (1 /cos(t)*log(u));
end
end

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

采纳的回答

Andrew Newell
Andrew Newell 2011-3-22
Here is some code that solves the problem and also illustrates an important programming principle: when you create a function, input the variables that you want to experiment with.
Here is a function that outputs a nonhomogeneous random sequence of times.
function S = nonhomogeneous(lambda0,lambda,T)
t = 0;
I = 0;
S = [];
u = rand;
t = t - log(u)/lambda0;
while t <= T
u = rand;
if (u <= lambda(t)/lambda0)
I = I+1;
S(I) = t;
end
u = rand;
t = t - log(u)/lambda0;
end
Save this in a file nonhomogenous.m. Here is a script to run it:
lambda0 = 50; % Needs to be much larger than T to get several numbers
T = 1;
lambda = @(x) lambda0*cos(x); % This is lambda0(t)/lambda0
S = nonhomogeneous(lambda0,lambda,T);
subplot 121
plot(S,lambda(S),'.')
xlabel('t')
ylabel('lambda(t)')
lambda = @(x) lambda0*sin(x);
S = nonhomogeneous(lambda0,lambda,T);
subplot 122
plot(S,lambda(S),'.')
xlabel('t')
I'm not sure what you mean by "spikes". Maybe you really want to plot something like
hist(S)
  3 个评论
显示 1更早的评论
Andrew Newell
Andrew Newell 2011-3-22
f(t) is lambda(t)/lambda. I should have written it this way for clarity, so I have edited the above code. Note that lambda0 is the maximum value for lambda in the code.
Susan
Susan 2011-3-22
Thanks very much :)

请先登录，再进行评论。

更多回答（2 个）

Andrew Newell
Andrew Newell 2011-3-22
@Susan, thank you for formatting the code nicely. I wish more people would in their questions!
I notice your function doesn't have any outputs and your inputs are redundant. If you change the top line to
function [s,i] = nonhomogeneous
do you get what you want?
EDIT: I see several problems:
t = 010; % t=10
instead of
t = 0;
Also,
T = length(cos(t));
gives T=1 (is that what you're after?)
s = [0,t];
initializes your array s to [0 10] instead of [].
Yet another: there is no lambda. Instead, you're using
max(cos(t))
which is only the maximum for the current values of t, not 1 as it should be.
  1 个评论
Susan
Susan 2011-3-22
Thanks very much Ive changed it accordingly but still I dont get what I am expecting.. I am meant to make the lambda function e.g sin or cos and then to have a graph at the end with that function and the higher the curve of the function the more spikes and vice versa..
the plot is not working at all.

请先登录，再进行评论。

Jan Jensen
Jan Jensen 2011-3-22
I've read the presentation, and it was quite interesting. Not sure I understand all the implications, though.
Please post your definition of the lambda function as well as the code.
A good starting point could also be to implement the simulation of Poisson process with constant decay (lambda = constant). The one on slide 5. The dimensions of the output variable S should be the same. If the plot of S for lambda=constant makes sense you have made a great step forward in debugging your code, including your plotting commands :)
  1 个评论
Susan
Susan 2011-3-22
Hey Jan..
I tried it with constant and it is workingm I just changed few things but its working finally :)

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Creating and Concatenating Matrices 的更多信息

标签

产品

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by