STA510 Statistical modelling and simulation, autumn 2020.

Mandatory assignment 1

Deadline: Friday September 19th at 12:00 (noon, Norwegian time).

Read carefully through the information about the page on Mandatory assignments on Canvas.

Notice in particular that the assignments should be solved individually.

Hand in on Canvas. Submit two files, one pdf-file with a report containing the answers to the

theory questions, and one file including the R-code. The first line of the R-code file should be:

rm(list=ls()) . Check that the R-code file runs before you submit it. Use comments in the

R-code to clearly identify which question each part of the R-code belong to. Also try to add some

comments to explain important parts of the code. The file ending of the R-code file should be .R

or .r. The report can be handwritten and scanned to pdf-file, or written in your choice of text

editor and converted to pdf. Cite the sources you use.

Problems marked with an R should be solved in R, the others are theory questions.

Problem 1:

In this problem, we consider a random variable X with probability density function given by

f(x) = (12 ?x4, if ? 1 ≤ x ≤ 10, otherwise

(1)

a) Verify that this is a proper probability density function.

Find the mean E(X) and the variance V ar(x).

Also find E.

b) Find the cumulative distribution function F(x) associated with X, and use this to calculate

P(X < 0) and P(?0.5 < X < 0.5).

Find also the inverse cumulative distribution function F

?1

(u) (Hint: You should get a

quadratic equation which generally have two solutions, but only one of these solutions

(which?) gives sensical results here.)

c)R Write a function (with single argument n) which produces n random numbers with density

(1) using the inverse transform method. Ideally, the function should not contain any forloops.

Check that your function produces correct results by comparing a histogram (use argument

probability=TRUE) of many random variables to the density (1).

d) Using a uniform(-1,1) proposal distribution g, do the calculations required to obtain an

accept-reject method for drawing random numbers with density (1).

1

e)R Write a function (with single argument n) which produces n random numbers with density

(1) using the accept-reject method.

Check the results in the same way as in point c.

Which method would you prefer in this case - accept-reject or inverse transform? The

function system.time is useful in this regard.

f)R Check the answers you got for expectations, variance and probabilities in points a) and b)

above by simulation. For E(X), also find a simulation-based 95% confidence interval. Use

either the inverse transform or accept-reject algorithm for generating random numbers with

density (1).

Problem 2:

In this problem, we let X1, X2, . . . , Xn be iid N(μ, σ2

) random variables.

a) Which distribution does A =Pn

i=1 Xi have?

Which distribution does B =1nPn

i=1 Xi have?

Which distribution does C =√n(B?μ)σ

have?

b)R Check the results you got in a) by simulating many independent replications of each of

A, B, C. You may use e.g. n = 100, μ = 1, σ = 2.

Now, think of X1, X2, . . . , Xn as random outcomes of some stochastic simulation, and you interested

in knowing the mean μ which you estimate using B above. Through initial runs of the

simulation, you have worked out that σ = 2 (at least approximately, and we will use this number

in subsequent calculations). To get the desired precision in the estimate B, you choose n so that

a 95% confidence interval is no longer than 0.1.

c) How many simulations n do you need to obtain such a level of precision in the mean estimate?

d)R For the n obtained in c) and μ = 1, simulate many outcomes of X1, X2, . . . , Xn. Based on

these outcomes, compute approximate sample mean confidence intervals on the form

where S is the sample standard deviation (e.g. obtained using R function sd()). How often

does the interval contain the true mean? Redo this calculation with n = 10. What do you

find?

Problem 3:

In this last problem, we consider the Kissler et al (2020) SEIR model. Please read the documentation

(available on Canvas) of the model carefully, and download the appropriate R-file

seir model.R before starting to solve this exercise. Throughout this problem, we consider a

population of one million persons - be sure to scale results appropriately.

2

a)R Using the default parameters, i.e. by running the model using

> f <- seir sim()

find the following quantities from the output f:

-The number of persons in state S (susceptible) on March 15th 2021.

-The maximum number of persons in critical care (CC) over the period, and also on which

date this occurs.

-The average number of persons in critical care in year 2021.

b)R In this point, we turn of social distancing by using the argument upper thresh=10000, i.e.

> f2 <- seir sim(upper thresh=10000)

-Find the maximum number of persons in critical care over the period.

-Using the definition that an epidemic ends when fewer than 1 per 1 million persons are in

the E-state, how long does it take to reach the end of the epidemic? (be careful here, as

the E-state is also zero before the epidemic starts...)

-Which proportion of the population has been infected (i.e. in RR, RH or RC) when the

epidemic ends in this model (herd immunity)?

Now suppose you are uncertain of the reproduction number R0 and the influence of seasonality

on the infection dynamics. To incorporate this uncertainty, you model the relevant parameters as

follows:

? R0max ～ uniform(1.8, 3.0).

? season ampliude ～ uniform(0, 0.4).

R0max and season ampliude are assumed to be independent. The remaining parameters are set

equal to their default values.

c)R Simulate many (at least 100) SEIR models with different random values of R0max and

season ampliude drawn from the above distributions. Based on these simulations; find the

following:

-The mean and standard deviation, the 90%, 95% and 99% quantiles of the maximum

number of required critical care beds during 2020.

-Calculate the probability that more than 100 critical care beds per 1 million persons are

required.

-Make a graphical representation of the maximum number of required critical care beds

during 2020. Is this distribution well approximated by a normal distribution?

d)R -Make scatterplots of the maximum number of critical care beds required in 2020 against

R0max and season ampliude.

-Would you say there is a linear relation between the required critical care beds and the

two parameters.

-Which of these parameters are influencing the required number of critical care beds?

-Calculate the correlations between the required number of critical care beds and the two

parameters.

3

版权所有：留学生编程辅导网 2018 All Rights Reserved 联系方式：QQ:99515681 电子信箱：99515681@qq.com

免责声明：本站部分内容从网络整理而来，只供参考！如有版权问题可联系本站删除。