Simulating computation of the expected value of random variable Y

I have been given the following task: Compute the expected value of Y=e^{-X}(X is uniform between 0 and 1) with a simulation in R. Plot the expected value as a function of the number of simulations where n is an integer between 1 and 10000 The pdf of this function is: f(y) = 1/y, for 1/e < y < 1.

The formula of finding expect value is of course: E[Y] = integrate(y * 1/y dy)

How do you simulate something likes this ? I expect you draw random sample between (1/e < y < 1) but the pdf of the distribution has different probabilities depending on what you draw it seems.

I thought about using the "sample" or "runif" functions but i can't figure out how to make those functions work with different probabilities.

You could use the fact that the random variable U = F(Y) is uniform, where F is the cumulative density function of the random variable Y (with pdf 1/y). You then have that Y = F^-1(U). This means that you can sample from a uniform variable and then transform it through F^-1(U) to get a sample from Y. You can then take the mean of your sample. This is known as inverse sampling transformation.

For your example you have F(y) = ln(y) + 1 and F^-1(u) = exp(u - 1). It is then easy to get a sample:

n = 1000
u = runif(n)
y = exp(u - 1)
mean(y)
0.6342477

which is very near the true mean of 0.6321206 (1 - 1/e).

EDIT

To see how the estimated mean changes with how many samples you simulate you could do something like this:

sample_y = function(n){
  u = runif(n)
  y = exp(u - 1)
  mean(y)
}

n = seq(10, 20000, 10)
res = sapply(n, sample_y)
ts.plot(res)

It stabilizes around the true mean very quickly and the variation around the mean gets smaller and smaller as n grows.