\documentclass[11pt]{article}
\usepackage{color} %used for font color
\usepackage{amssymb} %maths
\usepackage{amsmath} %maths
\usepackage{epsfig}
\topmargin-.5in
\textwidth6.6in
\textheight9in
\oddsidemargin0in
\def\ds{\displaystyle}
\def\d{\partial}
\def\bH{{\bf H}}
\def\bM{{\bf M}}
\def\bB{{\bf B}}
\def\bm{{\bf m}}
\def\bS{{\bf S}}
\def\ve{{\varepsilon}}
%\def\fs{\scriptsize}
\def\fs{\footnotesize}
\def\reals{\mathbb{R}}
\def\complex{\mathbb{C}}
\def\scriptO{{{\it O}\kern -.42em {\it `}\kern + .20em}}
\begin{document}

\centerline{\large \bf Math 540: Project 1}

\bigskip
\noindent
Due Thursday, January 30

\bigskip
\noindent
1. Consider the covariance function $C(x,y) = \min(x,y)$ on the interval $[0,1]$ along with the mean function $\bar{\alpha}(x) = \frac{1}{8}(x+1)^2$.  Compute the first $N=3$ eigenvalues and eigenfunctions and plot the eigenfunctions.  Now compute and plot 1000 realizations of the Karhunen-Loeve expansion truncated at $N=3$ with random variables $Q_n \sim N(0,1)$.  You can compute the random variables using the command {\tt randn.m}. 

For your 1000 realizations, use {\tt histnorm.m} to plot a histogram scaled to unity of the realizations of $\beta(\bar{x},\omega) = \sum_{n=1}^3 \sqrt{\lambda_n} \phi_n(\bar{x})Q_n(\omega)$ at $\bar{x} = 0.5$.  On the same figure, plot a kernel density estimate (KDE) of the distribution and a normal approximation to it.  Is the distribution unbiased and consistent with your choice of $\sigma^2=1$? Finally, repeat this plot with $N=6$ to see if your expansion has converged.

\medskip
\noindent
Note: Use the property that 
$$
  \mbox{\rm var}\left(\sum_{i=1}^n a_i X_i \right) = \sum_{i=1}^n a_i^2 \mbox{\rm var}(X_i)
$$

\bigskip
\bigskip
\noindent
2. For the same mean function $\bar{\alpha}(x) = \frac{1}{8}(x+1)^2$, now consider the exponential covariance function $C(x,y) = e^{-|x-y|/2L}$ on the interval $[-1,1]$ with the correlation length $L=5$.  The eigenvalues are 
$$
  \lambda_n = \left\{ \begin{array}{lcl} \frac{2L}{1 + L^2 w_n^2} & , & \mbox{\rm if $n$ is even,} \\
                                  \noalign{\smallskip}  \frac{2L}{1 + L^2 v_n^2} & , & \mbox{\rm if $n$ is odd,} 
                         \end{array} \right.
$$
and the eigenfunctions are
$$
  \phi_n(x) = \left\{ \begin{array}{lcl} \frac{\sin(w_n x) }{\sqrt{ 1 - \frac{\sin(2w_n)}{2 w_n}}} & , & \mbox{\rm if $n$ is even,} \\
                                  \noalign{\medskip} \frac{\cos(v_n x) }{\sqrt{ 1 + \frac{\sin(2v_n)}{2 v_n}}} & , & \mbox{\rm if $n$ is odd.}
                                  \end{array} \right.
 $$
 Here $w_n$ and $v_n$ are solutions of the transcendental equations
 $$
   \left\{ \begin{array}{lcl} Lw + \tan(w) = 0 & , & \mbox{\rm for even $n$, } \\
   \noalign{\smallskip}
                                             1 - L v \tan(v) = 0 & , & \mbox{\rm for odd $n$.}
   \end{array} \right.
$$
Compute the first odd and even roots $v$ and $w$.  You can get good initial guesses by plotting the functions and zooming to approximate the roots.  You can obtain accurate values using the MATLAB command {\tt fzero.m}.  Compute the first two eigenvalues and eigenfunctions and plot the eigenfunctions.  Now compute and plot 1000 realizations of the Karhunen-Loeve expansion truncated at $N=2$ with random variables $Q_n \sim N(0,1)$. 


\bigskip
\bigskip
\noindent
3. Consider again the covariance function $C(x,y) = \min(x,y)$ on the interval $[0,1]$ and a Karhunen-Loeve expansion truncated at $N=3$.  You can ignore the contributions from the random variables $Q_n$.  Your objective is to  approximate the mean $\bar{\alpha}(x)$ using the representation $\bar{\alpha}(x,q) = q_0 + q_1 x + q_2 x^2$, where $q = [q_0, q_1, q_2]$ are parameters to be estimated using the data in  \verb#KL_data.txt#.   In the file, the first column is the $x$-values at the 21 points $x_j = 0.05( j-1)$ and the second column is values of 
$$
\begin{array}{rl}
   y_j & = {\ds \alpha(x_j,q) + \varepsilon_j} \\
          \noalign{\smallskip}
         & = {\ds \bar{\alpha}(x_j,q) + \sum_{n=1}^3 \sqrt{\lambda_n} \phi_n(x_j) + \varepsilon_j }  
\end{array}
$$ 
where $\varepsilon_j$ is added measurement noise.  

You should use optimization software to solve the minimization problem
$$
  q_{opt} = \mbox{\rm arg}\min_{q} \sum_{j=1}^{21} \left[ y_j - \alpha(x_j,q) \right]^2.
$$
Plot your approximation to $\alpha(x,q)$ along with the data. 

You can see an example regarding the use of {\tt fminsearch.m} in Example 7.16 at the website \verb#https://rsmith.math.ncsu.edu/UQ_TIA/CHAPTER7/index_chapter7.html#.  We will revisit this problem again, to estimate the variances $\sigma_n$ for the random variables $Q_n$, once we have discussed frequentist and Bayesian inference.


\end{document}

