First actuarial example

I found this little example in an actuarial book using R.

set.seed(1)
X <- rexp(200, rate=1/100)
print(X[1:5])

I know that rexp generates random deviates.

How can I replicate this example in Julia? Using Distributions.jl but then what?

I think you are looking for

julia> using Distributions

julia> rand(Exponential(100), 100)

although double check the second parameter - I just guessed that it is the inverse of what R is using based on the name of the kwarg in R and the docstring for Exponential:

help?> Exponential

  Exponential(θ)

  The Exponential distribution with scale parameter θ has probability density function

  f(x; \theta) = \frac{1}{\theta} e^{-\frac{x}{\theta}}, \quad x > 0

  Exponential()      # Exponential distribution with unit scale, i.e. Exponential(1)
  Exponential(θ)     # Exponential distribution with scale θ
  
  params(d)          # Get the parameters, i.e. (θ,)
  scale(d)           # Get the scale parameter, i.e. θ
  rate(d)            # Get the rate parameter, i.e. 1 / θ
2 Likes

I tried with:

using Distributions
Random.seed!(1)
X = rand(Exponential(200), 100)
X[1:5]

but the results are different in R:
[1] 75.51818 118.16428 14.57067 13.97953 43.60686

Random.numbers are mostly just that - random. You can create a bit of consistentcy by setting a seed like you did but in general the random numbers will be implementation-specific. How Julia and R do random numbers is likely different and therefore you are seeing different output despite setting the same seed.

You could:

  • generate and store a series of random numbers from R and then use that wherever you need it
  • actually call R from Julia with RCall

P.S. welcome to Julia! Do check out JuliaActuary.org which has some examples and tutorials too.

2 Likes

Sorry you misread my post, but I made it easy for you to do so.

rand(Exponential(theta), n)

Draws n random numbers from an Exponential distribution with scale parameter theta. I interpreted Rs rate kwarg to be the inverse of the scale parameter in Distributions.

So when you do Exponential(200) you are likely sampling from a distribution that’s the equivalent of rate=1/200 in R. You should do rand(Exponential(100), 200)

(this is independent to Alex’s good point that random number generators are generally not comparably across programming languages so you shouldn’t expect to get the exact same numbers back from Julia and R even with the same random seed)

1 Like