M=c(1,2,3)
Ti=10  #number of times run simulation
c=0.2
A=5
d=10
d0=5
s=0.5
e=0.1
p=0.6
ER=e/A
C1 = c*A*exp(-d/d0) #1 to 3
C2 = c*A*exp(-d*s/d0) #1 to 2
C3 = c*p*A*exp(-(1-s)*(d/d0)) # 2 to 3
j=c(0:1)
k=c(0:1)
l=c(0:1)
prb=prb[j]*prb[k]*prb[l]

for(i in 1:M[3]){         #Location

  for(k in 1:Ti){          #Time length

M[3]=prb                    #inital colonization probability
if (M[3] > 0) {
  prb[j]=M[3]*ER           #then established
} else if (C1 > 0) {      #if not colonized/established
  prb[k]=C1               #prob 1 to 3
} else if (C3 > 0) {
prb[l]=C3                 #prob 2 to 3
}  else (C1+C3+M[3] < 0)
    (prb=0)
  }
}
prb