Hello, dear friend, you can consult us at any time if you have any questions, add WeChat: daixieit


STAT 440 Spring 2022 Homework 2

2022

 

This assignment is due on February 18, 2022 at 11:59 PM EST. To receive credit for an answer, you must concisely show your work and/or explain your reasoning. Please include all computer codes that support your work in the Canvas submission. Please follow the homework policy stated in the syllabus.

Problems marked (A) require you to derive analytic mathematical expressions for the solution, whereas problems marked (C) require you to write computer code. For each problem marked (C) , if you are writing your solution in R, please set your RNG seed using set.seed(440) at the beginning of your code solution.

 

1. A random variable X ∈ (0, ∞) has Gamma distribution, X ∼ Gamma(α , β), with

f (x) = xα − 1e −βx .

In this problem, you will investigate exponential tilting, a special way to choose proposal distributions for rejection and importance sampling. Let θ ∈ (0, β) and α > 1. The exponentially tilted distribution f θ(x) has the form:

f θ(x) = f (x).

Let Yθ be a random variable with the density fθ above.

Derivations:

E[eθX] = 1 −  −α .

ii. (A) Show that Yθ ∼ Gamma(α , β − θ).

iii. (A) Recall that E[X] = α/β. Let y∗ > α/β. Find the value of θ such that E[Yθ] = y∗ .

(b) Suppose we wanted to calculate P(X ∈ [a, b]) for some constants a, b ∈ R where a < b:

i. (A) Write the expression for the Monte Carlo estimate of P(X ∈ [a, b]) using iid samples from Gamma(α , β).

ii. (A) Write the expression for the Monte Carlo estimate of P(X ∈ [a, b]) using importance sampling with iid samples from Gamma(α , β − θ), where θ is chosen such that: E[Yθ] = (a+b)/2

iii. (A) Suppose a ≫ α/β. Which method do you think will perform better and why?

2. (Note: this problem is an application of the method you derived in problem 1; we recommend doing that problem first.)

In climate science, rainfall is frequently modeled using Gamma distributions. In Figure 1, we have the distribution of rainfall in mm from multiple geographic locations on the same day in Arkansas and Oklahoma. For more details, see Amburn, S. A., Lang, A. S. I. D., and Buonaiuto, M. A. (2015), “Precipita- tion Forecasting with Gamma Distribution Models for Gridded Precipitation Events in Eastern Oklahoma and Northwestern Arkansas”, Weather and Fore-

casting, 30(2), 349-367 (https://journals.ametsoc.org/view/journals/ wefo/30/2/waf-d-14-00054_1.xml). We see that the rainfall is reasonably ap- proximated with a Gamma(4.26, 6.07) model:

 

 

Figure 1: Distribution of Daily Rainfall Sampled from Arkansas and Oklahoma (Amburn et al, 2015)

Our goal is to estimate the probability that the rainfall in any sampled area is between 2mm and 3mm. Let X ∼ Gamma(4.26, 6.07)

(a) (C) Calculate the empirical (i.e. observed) frequency of rainfalls between

2mm and 3mm in the data.

(b) (C) Calculate the exact probability P(X ∈ [2, 3]) under the Gamma model using the pgamma function in R (or equivalent function in another language of your choice).

(c) (C) Estimate P(X ∈ [2, 3]) using the method in 1.b.i with n = 1000. Cal- culate the standard error for this method using the central limit theorem approximation.

(d) (C) Estimate P(X ∈ [2, 3]) using the method in 1.b.ii with n = 1000.            (e) (C) Of the two different methods, which method would you prefer, and why?

3. A Poisson process is a counting process, a special collection of random variables {N(t)}t≥0 that’s a function of time t ∈ [0, ∞) and a constant rate λ ∈ (0, ∞). The process models events that occur at random times and counts the number of events N(t) that have occurred up to time t. The Poisson process at time t has the following mass function for non-negative integers n:

P(N(t) = n) =  exp ( −λt)

In this problem, you can use the following properties about Poisson processes without proof:

• N(0) = 0 and N(t) is strictly increasing as t increases.

• For any two time points a < b:

N(b) − N(a) ∼ Poisson(λ(b − a))

• The Poisson process has independent increments, i.e. if a < b ≤ c < d, then: N(b) − N(a) ⊥ N(d) − N(c)

• If N(t) = n, then the times at which the n events occur in the interval [0, t] are iid Unif(0, t).

• If two sequential count events occur at times t1 and t2, then t2 − t1 ∼ Exp(λ). Suppose we observe the Poisson process over a time interval [0, T].

(a) (A) Show that:

N(t) | N(T) = N ∼ Binomial N, 

(b) (A) Using the previous result, show that N(T) is a sufficient statistic for λ .

(c) (A) Show that the estimator:

 ≜ N(T)

is an unbiased estimator for λ. Would Rao-Blackwellization improve this estimator? Why or why not?

4. (Note: this problem is an application of the method you derived in problem 3; we recommend doing that problem first.)

The data set arrival_times.csv lists arrival times for students to a testing center. All times are in hours, where t ∈ [0, 8] corresponds to times during a normal 9-to-5 day. We plot the data in Figure2. In this problem, we’ll model these times using a Poisson process.

 

Figure 2: Students arriving at the testing center over time

(a) (C) Estimate  using the method in part 3.b.

(b) (C) Write a function to simulate from the Poisson process with parameter  ,

and plot one realization of this Poisson process.


(c) (C) Perform the following hypothesis test using MC hypothesis testing: H0 : λ = 5 students per hour,   H1 : λ > 5 students per hour

Calculate the p-value using 1000 MC samples, and interpret the p-value in context.

(d) (C) Let γ be the average time between when students arrive at the testing center. Perform the following hypothesis test using MC hypothesis testing:

H0 : γ = .2 hours,   H1 : γ > .2 hours

Calculate the p-value using 1000 MC samples, and interpret the p-value in context.

5. We want to compute E[exp{−X2}] where the random variable X follows a non- central t distribution T (ν , µ , σ2). This distribution can be represented as

ε     X = µ +σ

ε ∼ N(0, 1),    ξ ∼ χ ,

where ε and ξ are independent. We can simulate random numbers from T (ν , µ , σ2):

• Step 1: Simulate Y = ξ/ν, where ξ ∼ χ .

• Step 2: Given Y from Step 1, simulate X = µ +σε/ where ε ∼ N(0, 1).

By running this simple algorithm n times we obtain an iid sample (Xi , Yi), i = 1,..., n, allowing us to estimate E[exp{−X2}] in two different ways. The first estimator is the empirical average (i.e., the naive Monte Carlo approach):

1  n

 

The second estimator is the Rao–Blackwellized average:

1  n

 

(a) (A) Derive the analytic form of Rn. You can find the final answer in our

lecture notes, but you need to show the derivation to get full credit here.

(b) (C) Compare the empirical distribution and variance of En  with Rn  for (n, ν , µ , σ) = (1000, 6, 2, 0.4) based on 10,000 independent replications. Note that the codes in our lecture notes are based on only one replication.