统计代写 | STAT 5361 Statistical Computing Exercises

3.5.1 Cauchy with unknown location.

Consider estimating the location parameter of a Cauchy distribution with a known scale parameter. The density function is$$\begin{array}{c}f(x;\theta )=\frac{1}{\pi [1+(x-\theta {)}^2]},x\in R,\theta \in R.\end{array}$$

Let $X_1,\dots ,X_n$

 be a random sample of size $n$

 and $\ell \left(\theta \right)$

 the log-likelihood function of $\theta$

 based on the sample. – Show that$$\begin{array}{cc}\ell \left(\theta \right)& =-n\ln \pi -\sum _{i=1}^n\ln [1+(\theta -X_i{)}^2],\\ {\ell }^{\prime }\left(\theta \right)& =-2\sum _{i=1}^n\frac{\theta -X_i}{1+(\theta -X_i{)}^2},\\ {\ell }^{\prime \prime }\left(\theta \right)& =-2\sum _{i=1}^n\frac{1-(\theta -X_i{)}^2}{[1+(\theta -X_i{)}^2{]}^2},\\ I_n\left(\theta \right)& =\frac{4n}{\pi }{\int }_{-\infty }^{\infty }\frac{x^{2}dx}{(1+x^2{)}^3}=n/2,\end{array}$$

where $I_n$

 is the Fisher information of this sample. – Set the random seed as $20180909$

 and generate a random sample of size $n=10$

 with $\theta =5$

. Implement a loglikelihood function and plot against $\theta$

. – Find the MLE of $\theta$

 using the Newton–Raphson method with initial values on a grid starting from $-10$

 to $30$

 with increment $0.5$

. Summarize the results. – Improved the Newton–Raphson method by halving the steps if the likelihood is not improved. – Apply fixed-point iterations using $G\left(\theta \right)=\alpha {\ell }^{\prime }\left(\theta \right)+\theta$

, with scaling choices of $\alpha \in \{1,0.64,0.25\}$

 and the same initial values as above. – First use Fisher scoring to find the MLE for $\theta$

, then refine the estimate by running Newton-Raphson method. Try the same starting points as above. – Comment on the results from different methods (speed, stability, etc.).

3.5.2 Many local maxima

Consider the probability density function with parameter $\theta$

:$$\begin{array}{c}f(x;\theta )=\frac{1-\cos (x-\theta )}{2\pi },0\le x\le 2\pi ,\theta \in (-\pi ,\pi ).\end{array}$$

A random sample from the distribution is

x <- c(3.91, 4.85, 2.28, 4.06, 3.70, 4.04, 5.46, 3.53, 2.28, 1.96,
       2.53, 3.88, 2.22, 3.47, 4.82, 2.46, 2.99, 2.54, 0.52)

• Find the the log-likelihood function of $\theta$
• Find the method-of-moments estimator of $\theta$
• Find the MLE for $\theta$
• What solutions do you find when you start at ${\theta }_0=-2.7$
• Repeat the above using 200 equally spaced starting values between $-\pi$

3.5.3 Modeling beetle data

The counts of a floor beetle at various time points (in days) are given in a dataset.

beetles <- data.frame(
    days    = c(0,  8,  28,  41,  63,  69,   97, 117,  135,  154),
    beetles = c(2, 47, 192, 256, 768, 896, 1120, 896, 1184, 1024))

A simple model for population growth is the logistic model given by$$\frac{dN}{dt}=rN(1-\frac{N}{K}),$$

where $N$

 is the population size, $t$

 is time, $r$

 is an unknown growth rate parameter, and $K$

 is an unknown parameter that represents the population carrying capacity of the environment. The solution to the differential equation is given by$$N_t=f\left(t\right)=\frac{K{N}_0}{N_0+(K-N_0)\exp (-rt)},$$

where $N_t$

 denotes the population size at time $t$

.

• Fit the population growth model to the beetles data using the Gauss-Newton approach, to minimize the sum of squared errors between model predictions and observed counts.
• Show the contour plot of the sum of squared errors.
• In many population modeling application, an assumption of lognormality is adopted. That is , we assume that $\log N_t$