\documentclass[11pt]{article}

\usepackage{color}
\usepackage{graphicx}
\usepackage{latexsym}
\usepackage{wrapfig}
\usepackage{url}

\setlength{\topmargin}{0in}
\setlength{\headheight}{0in}
\setlength{\headsep}{0in}
\setlength{\textheight}{9.0in}
\setlength{\textwidth}{6.5in}
\setlength{\oddsidemargin}{0in}
\setlength{\evensidemargin}{0in}
\setlength{\parindent}{0.0in}
\setlength{\parskip}{0.1in}
\renewcommand{\baselinestretch}{1.2}


\begin{document}
\vspace*{-8mm} 
{\Large\centerline{Reference Manual for Package ``sparse-msrf''}}
\vspace{3mm}
{\large{\centerline{J. Ray and J. Lee}
\centerline{Sandia National Laboratories, Livermore, CA }
\centerline{\small \{jairay,jlee3\}@sandia.gov}
}
\vspace{2mm} 

%%%%
%% SAND2014-16638 O
%%%%

\section{Introduction}
This is a manual for the package ``sparse-msrf'', a Matlab package for
estimating fossil-fuel $\rm CO_2$ (ffCO2) emissions from sparse measurements
of its concentration at a limited number of locations. This software
was used in the study documented in~\cite{sandrept,
paper:msrf,paper:theory}. It is suggested that you read either~\cite{sandrept}
or \cite{paper:msrf,paper:theory} before delving into this software
package. 

Estimation of fields (defined on a grid), from sparse data, cannot be
performed naively i.e. by assuming that each grid-cell is an
independent parameter to be estimated - there is simply not enough
observations to estimate all the cells. If the field is smooth, and
adjoining grid-cells are constrained to assume field values which are
not too different, one can use Gaussian Random Field~\cite{hd061a} to
model them and drastically reduce their dimensionality (making it
feasible to estimate the lower dimensional model from limited
data). When fields are not smooth, one can model them as a (linear)
combination of proxies that contain similar patterns or one can use
basis sets with compact support e.g., wavelets. A wavelet
encoding/model of a field is not low-dimensional, unless one can
``mask'' (or eliminate) a large fraction of the wavelets using a
guess/easily observed proxy of the field being modeled (which then
results in a random field model). Since wavelets are a
multi-resolution basis, one can called a masked collection of wavelets
a Multiscale Random Field model (MsRF).

Given the paucity of data, even the model may not be sufficiently low
dimensional. In such a case, one takes recourse to sparse
reconstruction (shrinkage regression) which identifies which wavelet
weights can be estimated from observations while turning off the
rest. This is most efficient when the parameters being estimated
(wavelet weights) are independent (alternatively, when the bases are
orthogonal); since wavelets are an orthogonal basis set, this hold true
for us.

The software contains the following capabilities:
\begin{itemize}
  \item The ability to search through a number of wavelet families to
  determine the sparsity of encoding a given non-stationary field

  \item The functionality to create a MsRF model by masking specific
  wavelets in a wavelet encoding of a field. The field being used to
  create the masking could be an emission field or an easily-observed
  proxy (as was in our case).

  \item Functionality to construct various random sampling matrices
  and compare their sampling efficiency versus the atmospheric
  transport model.

  \item A sparse reconstruction functionality, based on
  StOMP~\cite{06dt4a}, to estimate an emission field in an irregular
  (non-rectangular) domain while ensuring non-negativity (ffCO2
  emissions cannot be negative).
\end{itemize}

The software package can be easily adapted to conduct network design
or Observation System Simulation Experiments for ffCO2 estimation. It
can be modified to also estimate emission patterns for gases which can
be approximated as non-reactive tracers e.g. methane emissions from
landfills, when the number, and distribution of landfills is
ambiguous~\cite{13hj5a}.

\emph{Citation:} The proper way to cite this package is:
``J. Ray and J. Lee, ``sparse-msrf: A package for sparse modeling and estimation of
fossil-fuel CO2 emission fields'', \url{http://www.sandia.gov/~jairay/kfcs.html}.

\section{Description of the package}

We envision that users will use this package to (1) experiment with
wavelets to model non-stationary field (not just ffCO2 emission
fields) and (2) use the software to perform source or emission field
inversions for passive tracers. To that end, we also provide a
realistic transport model.

The contents of the package are separated into different directories.
\begin{itemize}
  \item \texttt{data/}: This directory contains miscellaneous data,
  required as input for the study in our paper.

  \begin{itemize}

    \item \texttt{vulcan.8day.1/}: Folder containing ffCO2 emissions
    from the Vulcan inventory for US~\cite{web:vulcan}, aggregated to
    1 degree and 8 day resolution. This is the ``ground truth''
    emission field and is used to select wavelet family for the MsRF.

    \item \texttt{nlights-0.1x0.1/}: A folder containing radiances of
    lights at night over North America. This is an easily-observed
    proxy for ffCO2 emissions and is used to choose which wavelets to
    mask.

    \item \texttt{edgar2005.1x1/}: Emission field from the EDGAR
    inventory~\cite{web:edgar} which is used to provide a
    constant-in-time approximation of time-variant ffCO2 emission
    fields to be estimated.

    \item \texttt{tower.obs}: Names of towers where ffCO2
    concentrations are measured.

  \end{itemize}

  \item \texttt{doc/}: Contains this manual and bibliography.

  \item \texttt{examples/:} This folder contains 5 examples that drive
  the modeling and algorithmic machinery of the package. 3 examples in
  (\texttt{examples/ex01/}) concern wavelet modeling of ffCO2 fields
  using a variety of wavelets, the sparsity of representation and
  generally exercise the machinery for representing a masked
  collection of wavelets. \texttt{examples/ex02} contains a large
  example that uses Haar wavelet and nightlight radiances to construct
  a MsRF, uses EDGAR to create a prior emission field and uses StOMP
  to perform the estimation. This code is quite complicated and ought
  to be learnt piecemeal. \texttt{examples/ex03/} contains code that
  compares the sampling efficiency of the transport model versus
  traditional random sampling matrices.

  \item \texttt{H-matrices/}: ffCO2 concentration at a space-time
  location $i$, ${\mathbf y}( x_i, t) = {\mathbf H}_i {\mathbf z}$
  where ${\mathbf z}$ is the emission field and ${\mathbf H}$ is the
  linear operator that serves as the transport model over North
  America, with emission at 1 degree and 8-day resolution. This
  directory contains the ${\mathbf H}_i$ for 30+ observation towers.

  \item \texttt{incoherence/}:  Folder with software machinery to
  compute random sampling matrices and their coherence with the
  wavelet bases constituting the MsRF.

  \item \texttt{msrf-selection/}: This folder contains the software
  machinery for representing a fields as a masked collection of
  wavelets of a number of families and orders. The masking field (in
  our case, a nightlight matrix) is an argument to the Matlab
  functions that perform the masking.

  \item \texttt{sparse-inversion/} : This folder contains the bulk of
  the code. It has functions that select the wavelets for masking
  using nightlights, use EDGAR to make the guess emission field, and
  perform estimation of the non-negative field (Steps I and II, as
  described in~\cite{paper:theory}, are clearly marked). All the
  functions are documented, but it is expected that you might need
  help. Mail Jaideep Ray.

  Do not try to understand the code by looking through
  \texttt{msrf-selection/} and \texttt{sparse-inversion/}; rather
  start with the examples and work your way down.
\end{itemize}

\section{Installation}

The package depends on two open-source Matlab toolboxes,
SparseLab~\cite{web:sparselab} and WaveLab~\cite{web:wavelab}. Install
them first. They contain their own installation manuals.

Uncompress the \texttt{sparse-msrf.tgz} file 

In each of the examples, you will have to set the path to the two
toolboxes before you can run them. These specifications are at the top
of the example Matlab scripts.

\section{Description of examples}

\texttt{examples/ex01/Fig02.m}: This example shows how different types
of wavelets provide wavelet coefficients of different average values
(and their variance) at different spatial scales, when encoding
annually-averaged ffCO2 emissions from Vulcan. Run this once so that
all useful function in the toolboxes are in your path.

Use this example to explore the encoding different wavelet families
and different orders inside the same family can provide. If you are
trying to model a different field e.g., methane emissions, chances are
that Haar wavelets (which we found to be best) is the wrong choice. Type
\texttt{help MakeONFilter} to determine the wavelet families and
corresponding order that are implemented in WaveLab. All user-settable
parameters are at the top of the Matlab script.

\texttt{examples/ex01/Fig23.m}: In any wavelet encoding, the bulk of
the wavelets are at the finest scales; the ability to mask off a large
fraction of them results in an enormous reduction of dimensionality (a
\emph{sparse} representation). This examples show how certain wavelet
families do not manage to yield a sparse representation; wavelet
weights at fine scale are all rather small, but not zero, and all
about the same magnitude, making it impossible to construct as masking
threshold.

Use this example to explore how the masking metric becomes absolutely
impossible as one progresses to smoother wavelets with larger support,
while the ffCO2 field has sharp discontinuities at the 1 degree
resolution scale.

\texttt{examples/ex01/Fig24.m}: The more wavelets we remove, the more
approximate the model and bigger the error when we try to model Vulcan
emissions with them. This code allows you to see how model fidelity
drops as one strives for a sparse model.


\texttt{examples/ex02/driver01.m}: This is a complicated example that,
as distributed, performs an inversion over a 3-Week period (Week 31 to
33, where a Week is 8 days). One needs a buffer of 2 Weeks ahead of the
period of interest. This is the place to start if you wish to delve
into the workings of the modeling and inversion code. The example will
take about 10 minutes to run.

Use this example to explore how to do a year-long inversion (see line
55). It'll take 4 hours and require 20 GB RAM; memory requirements
increase with the duration over which you are estimating the
emissions.  Also check how the inversion worsens when we do a bad job
of enforcing zero emissions outside US (reduce
\texttt{info.nCompSamples} on line 49)

The example uses StOMP, as implemented in
\texttt{SparseLab2.1-Core/Solvers/SolveStOMP.m} (in the SparseLab toolbox that
you install). The StOMP solver has some hardcoded parameters. Change
\begin{verbatim}
{[x_I, flag] = lsqr(@lsqrMat,y,OptTol,20);}
\end{verbatim}
on line 122 to
\begin{verbatim}
{[x_I, flag] = lsqr(@lsqrMat,y,OptTol,500);} 
\end{verbatim}
and recheck.


\texttt{examples/ex03/Fig21.m} contains an example where the
incoherence of the ${\mathbf H}$ matrix, with respect to the MsRF, is
compared to the incoherence achieved with standard random matrices
(which are generated as part of the code). The example is distributed
with incoherence estimation over a 2 Week period (to save on memory).

Explore how the incoherence changes throughout the year by testing for
bigger durations (not just 2 Weeks) and at different points in the
year (not just Week 21-23).

\section{Disclaimer of warranty}
We make no warranties, explicit or implicit, that the programs
contained in this collection are free of error, or are consistent with
any particular standard of accuracy, or that they will meet your
requirements for any particular application. They should not be relied
on for any purpose where incorrect results could result in loss of
property or personal injury. If you do use programs for any such
purpose, it is at your own risk. The authors disclaim all liability of
any kind, either direct or consequential, resulting from your use of
these programs.

The software is being released under a BSD license. Read LICENSE.TXT to 
understand the restrictions on the use and distribution of this software 
package.

\section*{Acknowledgments}
This work was funded by the Laboratory Directed Research \&
Development (LDRD) program at Sandia National Laboratories. 

Sandia National Laboratories is a multi-program laboratory managed and
operated by Sandia Corporation, a wholly owned subsidiary of Lockheed
Martin Corporation, for the U.S. Department of Energy's National
Nuclear Security Administration under contract DE-AC04-94AL85000.
\bibliographystyle{unsrt} \bibliography{manual}
\end{document}
