############################################
##  Proceso de Yule ########################
############################################

#Damos el punto inicial k>=1, el parámetro lambda=lam y el intervalo a simular [0,T]

simYule<-function(k,T,lam)
{
    Tsalto<-c()   #vector que guarda los tiempos de salto del proceso
    Nsaltos<-0
    ind=1          
    TAcum=0       #tiempo acumulado
    while(TAcum<=T)
    {
        Tsalto[ind]=TAcum+rexp(1,rate=lam*(k+Nsaltos))
        TAcum=Tsalto[ind]
        Nsaltos=Nsaltos+1
        ind=ind+1
    }
    
    
    plot(c(0,Tsalto),k:(k+Nsaltos),xlim=c(0,T),ylim=c(0,k+Nsaltos),type="s",xlab="Tiempo",
         ylab="Poblacion",main="Simulacion Proceso de Yule")
    return(Tsalto)
}

simYule(1,10,0.2) #Simulacion de Yule en [0,10], que inicia con 1 individuo y lambda=0.2
simYule(1,10,1) #Simulacion de Yule en [0,10], que inicia con 1 individuo y lambda=1
simYule(2,10,0.2) #Simulacion de Yule en [0,10], que inicia con 2 individuos y lambda=0.2


############################################
##  Proceso Explosivo ######################
############################################

Tsalto<-c()   
Nsaltos<-0
ind=1          
TAcum=0       
MaxS=1000000

while(Nsaltos<=MaxS)
{
  Tsalto[ind]=TAcum+rexp(1,rate=((Nsaltos+1)^2))
  TAcum=Tsalto[ind]
  Nsaltos=Nsaltos+1
  ind=ind+1
}

TAcum

plot(c(0,Tsalto),1:(Nsaltos+1),xlim=c(0,TAcum+0.1),ylim=c(0,10000),type="s",xlab="Tiempo",
     ylab="Poblacion",main="Simulacion Proceso Explosivo")