# Simulation the event times occuring in (0,T) of the delayed renewal process # with stationary increments such that T2 ~ U(0,1) T <- 100 t <- 0 I <- 0 # The jump times of the renewal process are stored in S S <- c() repeat{ u <- runif(1) # The distribution of T1 is given by 1 - sqrt(1 - T2) if(I == 0){ u <- 1 - sqrt(1 - u) } # A new interarrival time is added to time t t <- t + u if(t > T){ break } I <- I + 1; S <- cbind(S,t) } # Plot plot(c(0,S), seq(0, length(S), 1), type="s", xlab='t', ylab='N_t')