[R] Weibull point process
Torbjørn Ergon
t.h.ergon at bio.uio.no
Sat Sep 17 10:38:36 CEST 2011
On Thu, 15 Sep 2011 22:12:54 +1200, Rolf Turner
<rolf.turner at xtra.co.nz> wrote:
> On 15/09/11 19:24, Torbjørn Ergon wrote:
>> On Thu, 15 Sep 2011 09:47:35 +1200, Rolf Turner
>> <rolf.turner at xtra.co.nz> wrote:
>>> On 15/09/11 07:21, Torbjørn Ergon wrote:
>>>> Dear list,
>>>>
>>>> I'm looking for a function to generate (simulate) a random Weibull
>>>> point process. Can anyone help?
>>>>
>>>> Cheers,
>>>>
>>>> Torbjørn Ergon, University of Oslo
>>>
>>> Do you mean a renewal process with the inter-event times having
>>> a Weibull distribution? Should be trivial to code up, using
>>> rweibull.
>>>
>>> cheers,
>>>
>>> Rolf Turner
>>
>>
>> Thanks for your reply!
>>
>> Yes, in a way - but the parameters of the Weibull must change after
>> each event. I'm trying to simulate events that can happen multiple
>> times thru the life of an individual and where the hazard rate is age
>> dependent according to a Weibull hazard. Hence, after the first event
>> has occurred at time t1 the hazard rate should be h(x; t1,
>> shape,scale) = (shape/scale)*((x-t1)/scale)^(shape-1). This doesn't
>> seem trivial to do using 'rweibull' - but perhaps I'm missing
>> something (trivial).
>
> Unless I am terribly confused your hazard function is just the hazard
> function
> of a Weibull distribution ``starting off'' from t1 rather than
> starting off from 0.
> So a random variable having that hazard function would have the same
> distribution
> as t1 + X where X is Weibull with the given shape and scale.
>
> So you can just create the points of your process as the cumulative
> sum of
> a number of independent Weibull variates.
>
> Perhaps I'm not seeing things correctly. If so, would some wiser
> person please chip in and
> set me straight?
>
> cheers,
>
> Rolf Turner
Thanks again for your reply and sorry I was so bad at explaining this!
I'm trying to simulate a process where an event can happen multiple
times and where the hazard rate declines with age of the individual (not
just with time since last event) [I wasn't quite correct in saying that
"the parameters of the Weibull must change after each event", and it
should be x+t1 in the formula I gave].
After some thinking, I have come to that I (think I) can simulate this
by discretizing and drawing events from a Bernoulli with probablity
1-exp(-integrated hazard) for each time bin - see script below. An
alternative way of simulating this would be to draw X number of values
from a Weibull distribution where X is stochastic, but I don't know what
distribution X should have.
Does this make sense?
rweibull.point.proces = function(end.time, shape, scale, delta =
end.time/1000){
event.times = NULL
for(t in seq(0,end.time, by=delta)){
p = 1 - exp(-(((t+delta)/scale)^shape - (t/scale)^shape)) #
probability of an event occuring betweem time t and t+delta
z = rbinom(1,1,p)
if(z==1) event.times = c(event.times, t+delta/2)
}
return(event.times)
}
shape = 0.5
scale = 4
end = 20
weibull.hazard = function(x,k,L) (k/L)*(x/L)^(k-1)
par(mfrow=c(1,2))
x = seq(0,20,0.1)
plot(x,weibull.hazard(x,k=shape, L=scale), type="l")
y = rweibull.point.proces(end, shape, scale)
plot(y,y ,xlim = c(0,end), ylim= c(0,end))
Cheers,
Torbjørn Ergon
More information about the R-help
mailing list