\documentclass{article}

\usepackage{fullpage}
\usepackage{amsmath,amssymb}
\usepackage{bm}
\usepackage{dsfont}
\usepackage{enumitem}
\usepackage{xcolor}
\usepackage{url}
\usepackage{verbatim}
\usepackage{graphicx}
\usepackage{parskip}

\usepackage[colorlinks,linkcolor=red!80!black]{hyperref}
\usepackage[noabbrev,capitalize]{cleveref}

% Use one or the other of these for displaying code.
% NOTE: If you get
%  ! Package minted Error: You must invoke LaTeX with the -shell-escape flag.
% and don't want to use minted, just comment out the next line
\usepackage{minted} \BeforeBeginEnvironment{minted}{\begingroup\color{black}} \AfterEndEnvironment{minted}{\endgroup} \setminted{autogobble,breaklines,breakanywhere,linenos}

\usepackage{listings}


% Commands for questions and answers
\definecolor{question}{rgb}{0,0,1}
\newcommand{\ask}[1]{\textcolor{question}{#1}}
\newenvironment{asking}{\begingroup\color{question}}{\endgroup}

\crefname{section}{Question}{Questions}
\newlist{qlist}{enumerate}{1}
\setlist[qlist,1]{leftmargin=*, label=\textbf{[\thesection.\arabic*]}, ref={[\thesection.\arabic*]}}
\crefname{qlisti}{Question}{Questions}

\definecolor{answer}{rgb}{0,.5,0}
\newcommand{\ans}[1]{\par\textcolor{answer}{Answer: #1}}
\newenvironment{answer}{\par\begingroup\color{answer}Answer: }{\endgroup}

\newcommand{\red}[1]{\textcolor{red}{#1}}

\definecolor{points}{rgb}{0.6,0.3,0}
\newcommand{\pts}[1]{\textcolor{points}{[#1~points]}}
\newcommand{\onept}[1]{\textcolor{points}{[1~point]}}

\newcommand{\hint}[1]{\textcolor{black!60!white}{\emph{Hint: #1}}}
\newcommand{\meta}[1]{\textcolor{black!60!white}{\emph{#1}}}

\newcommand{\TODO}{\red{TODO}}

% misc shortcuts
\newcommand{\fig}[2]{\includegraphics[width=#1\textwidth]{#2}}
\newcommand{\centerfig}[2]{\begin{center}\includegraphics[width=#1\textwidth]{#2}\end{center}}

% Math
\DeclareMathOperator*{\argmax}{arg\,max}
\DeclareMathOperator*{\argmin}{arg\,min}
\DeclareMathOperator*{\E}{\mathbb{E}}
\newcommand{\bigO}{\mathcal{O}}
\newcommand{\bmu}{\bm{\mu}}
\newcommand{\bSigma}{\bm{\Sigma}}
\newcommand{\bX}{\mathbf{X}}
\newcommand{\by}{\mathbf{y}}
\newcommand{\eye}{\mathbf{I}}
\DeclareMathOperator{\indic}{\mathds{1}}
\DeclareMathOperator{\N}{\mathcal{N}}
\newcommand\R{\mathbb{R}}
\newcommand{\tp}{^\mathsf{T}}
\newcommand{\ud}{\,\mathrm{d}}
\DeclareMathOperator*{\Var}{Var}

\newcommand{\toth}[1]{^{(#1)}}
\newcommand{\ith}{\toth{i}}

\DeclareMathOperator{\Poisson}{Poisson}

%%begin novalidate  % overleaf code check gives false positives here
\makeatletter
% \abs{} uses auto-sized bars; \abs*{} uses normal-sized ones
\newcommand{\abs}{\@ifstar{\@Abs}{\@abs}}
\newcommand{\@abs}[1]{\left\lvert #1 \right\rvert}
\newcommand{\@Abs}[1]{\lvert #1 \rvert}

\newcommand{\norm}{\@ifstar{\@Norm}{\@norm}}
\newcommand{\@norm}[1]{\left\lVert #1 \right\rVert}
\newcommand{\@Norm}[1]{\lVert #1 \rVert}

\newcommand{\ceil}{\@ifstar{\@Ceil}{\@ceil}}
\newcommand{\@ceil}[1]{\left\lceil #1 \right\rceil}
\newcommand{\@Ceil}[1]{\lceil #1 \rceil}

\newcommand{\floor}{\@ifstar{\@Floor}{\@floor}}
\newcommand{\@floor}[1]{\left\lfloor #1 \right\rfloor}
\newcommand{\@Floor}[1]{\lfloor #1 \rfloor}

\newcommand{\inner}{\@ifstar{\@Inner}{\@inner}}
\newcommand{\@inner}[2]{\left\langle #1, #2 \right\rangle}
\newcommand{\@Inner}[2]{\langle #1, #2 \rangle}
\makeatother
%%end novalidate  % overleaf code check false-positives here


\begin{document}

\begin{center}
\Large
CPSC 440/550 Advanced Machine Learning (Jan-Apr 2024)\\
Bonus Assignment --
due Saturday April 27th at \textbf{11:59pm}
\end{center}


The assignment instructions are the same as for the previous assignments. If you do the assignment with a partner, please only hand in one copy for the group (using the appropriate Gradescope feature).






\section{Discrete Markov Chain Inference \pts{40}}

The function \texttt{grad\_chain()} in \texttt{main.py} makes a Markov chain model corresponding to the ``CS graduate Markov chain'' from class.

Note that, slightly contrary to class, this question will use the notation that the first step is index 0;
this will help avoid confusion between the math and zero-indexed Python code:
\[
    p(x_0,x_1,\dots,x_{d-1}) = p(x_0) \prod_{j=1}^{d-1}{p(x_j\mid x_{j-1})}
.\]

\begin{qlist}
\item \pts{6}
    Implement the function \texttt{MarkovChain.sample}, in \texttt{markov\_chain.py}, that uses ancestral sampling to sample sequences of a given length $d$. \ask{Hand in this code, and report the univariate marginal probabilities for time 30 using a Monte Carlo estimate based on 10,000 samples (from \texttt{main.py mc-sample}).}

\hint{\href{https://numpy.org/doc/stable/reference/random/generated/numpy.random.Generator.choice.html}{\texttt{rng.choice}} might be helpful.}

\begin{answer}\TODO\end{answer}

\item \pts{6}
    Write a function, \texttt{MarkovChain.marginals}, that uses the CK equations to compute the exact univariate marginals up to a given time $d$. Use \texttt{main.py mc-marginals} to find the exact marginals at time 30. \ask{Hand in this code, report the exact univariate marginals at time 30 (to three decimal places), and report how this differs from the marginals in the previous question.}
\begin{answer}\TODO\end{answer}

\item \pts{2}
    \ask{What is the state $c$ with highest marginal probability, $p(x_j = c)$, for each time $j$?}
\begin{answer}\TODO\end{answer}

\item \pts{9}
    Write a function, \texttt{MarkovChain.mode}, that uses the Viterbi decoding algorithm for Markov chains to find the optimal decoding up to a time $d$.
    \ask{Hand in this code and report the optimal decoding of the Markov chain up to time 50, and up to 100.}
\begin{answer}\TODO\end{answer}

\item \pts{2}
    \ask{Report all the univariate conditional probabilities at time 30 if the student starts in grad school, $p(x_{30} = c \mid x_0 = 2)$, for all $c$}.

\hint{you can do this by changing the input to the CK equations.}

\begin{answer}\TODO\end{answer}

\item \pts{4}
    \ask{Report for all $c$ the univariate conditional probabilities $p(x_4 = c \mid  x_{9} = 5)$ (``where you were likely to be 4 years after graduation, if you ended up in academia after 9 years'') obtained using a Monte Carlo estimate based on 10,000 samples and rejection sampling. Also report the number of samples accepted from the 10,000 samples.}
\begin{answer}\TODO\end{answer}

\item \pts{9}
    \ask{Give code implementing a dynamic programming approach to exactly compute $p(x_3 = c \mid  x_{6} = 5)$, and report the exact values for all $c$.}
\begin{answer}\TODO\end{answer}

\item \pts{2} \ask{Why is $p(x_j = \red{6} \mid x_{20} = \red{5})$ equal to zero for all $j$ less than $20$?}
\begin{answer}\TODO\end{answer}

\end{qlist}


\clearpage
\section{MCMC for Bayesian Logistic Regression \pts{20}}

\texttt{main.py mcmc-blogreg} loads a set of images of `2' and `3' digits.
It then runs the Metropolis MCMC algorithm to try to generate samples from the posterior over $w$, in a logistic regression model with a Gaussian prior.
Once the samples are generated, it makes a histogram of the samples for several of the variables.
The plots also show the MAP estimate, as a vertical line.\footnote{The division into ``positive,'' ``negative,'' and ``neutral'' variables is based on the sign of that MAP estimate.}

\begin{qlist}

\item \pts{4} \ask{Why might the samples coming from the Metropolis algorithm as run here not give a good approximation to the posterior?}
\begin{answer}\TODO\end{answer}

\item \pts{4} Modify the proposal used by the demo to $\hat{w} \sim \mathcal{N}(w,(1/10)^2 \eye)$ instead of $\hat{w} \sim \mathcal{N}(w, \eye)$. \ask{Hand in your modifications to the code and the updated histogram plot.}

\begin{answer}\TODO\end{answer}

\item \pts{4} Modify the proposal to use $\hat{w} \sim \mathcal{N}(w,(1/100)^2 \eye)$. \ask{Do you think this is probably sampling from the posterior distribution more or less effectively than the previous choice? (Briefly explain.)} (Don't base your answer on the test likelihoods, those are for the next part!)
\begin{answer}\TODO\end{answer}

\item \pts{8} The \texttt{BayesianLogisticRegression.different\_log\_likelihoods} method implements four different attempts at evaluating the predictive log-likelihood of this method on test data. Each one has a different interpretation. \ask{Explain each computation}; don't just directly translate the code into English, but say what it's computing in words. \ask{Which is the correct log-likelihood for the fully-Bayesian prediction?}

\end{qlist}


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\clearpage
\section{Graph Neural Networks \pts{40}}

While not nearly as common as tabular data, images, or natural language, another major modality of data that's seen a lot of interest particularly in the past five years or so is \emph{graph} data. There are multiple related problem setups, but here we wish to learn a \emph{function on graphs}: similar to how in image classification we might want to learn a function whose input is an image and output is a class label, here we want to learn a function whose input is the molecular graph of a peptide, and whose output gives some aggregate 3D properties of that molecule's structure, such as its length, ``sphericity,'' and inertia.

\begin{center}
\includegraphics[width=0.4\textwidth]{figs/peptides-struct-1.png}
~
\includegraphics[width=0.4\textwidth]{figs/peptides-struct-2.png}
\end{center}

This is the \texttt{Peptides-struct} dataset from the \href{https://arxiv.org/abs/2206.08164}{Long-Range Graph Benchmark}.
Each node represents a heavy atom in the peptide, and the (undirected) edges a molecular bond; nodes have pre-computed features associated with them.
(There are also edge features in this dataset, but we're ignoring them.)

We've pre-processed the dataset to include eigendecompositions of the graph Laplacian,
which will be useful for constructing Laplacian positional encodings in the \texttt{LapPENodeEncoder} class in \texttt{gnns.py}.
This gives some notion of ``location'' within the graph to each node,
similar to the trigonometric position features we discussed for sequence Transformers.
(We won't use this for the GCN, though you could, and it would help.)
You don't have to worry about this, but you can check out that class (and the code inside the \texttt{if False} block in \texttt{gnn\_utils.py}) if you're curious.

To work with graph data, you'll want to install the \href{https://pytorch-geometric.readthedocs.io/en/stable/}{PyTorch Geometric} library, with \texttt{pip install torch\_geometric}.
You won't have to worry about any internals here, but some of its helpers will handle some of the grunt work of working with this kind of data.
The first time you run \texttt{main.py gcn}, it'll download the dataset into the \texttt{data} folder; it's about 400 MB, so I didn't put it in the zip.

This dataset is bigger than the datasets we've used before (ie MNIST); it's still runnable on a decent laptop, but it might be a little annoying,
especially for the Transformer later.
You might prefer to use a GPU on Google Colab (which is free),
or some other machine with a CUDA-capable GPU.
To use Colab, go to \url{https://colab.research.google.com}, and request a T4 GPU under \texttt{Runtime $\to$ Change runtime type}.
Then run \texttt{\!pip install torch\_geometric},
open the Files tab on the left (the folder icon) and upload the contents of the code directory.
Then, after \texttt{import main},
you can use either \texttt{main.train\_gcn()} or \texttt{main.run("gcn")}.
Be careful with editing files on Colab, though; unless you save them to your Google Drive or similar,
they'll just vanish on you!
One thing you could do is write your code locally,
make sure it runs on a batch or two,
then upload to Colab for a final run to make sure it trains okay.

If you have CUDA set up for your GPU in pytorch, e.g. on Colab, it should use that automatically.
You can use \texttt{main.py use-cpu gcn}, \texttt{main.py use-mps gcn}, or \texttt{main.py use-cuda gcn} to force a device.
MPS is the GPU on recent Macs; on my machine, I get pretty mixed results with CPU usually faster on this workload, but your mileage may vary.
If you run into memory issues, try decreasing the batch size given to the data loader; if you have plenty of free memory, you could try increasing it to see if that's faster.

\texttt{main.py gcn} will run a Graph Convolutional Network on this data, using the \texttt{torch\_geometric} library's implementation.
This model is, while not state of the art or anything, ``pretty okay'' on this dataset.
This dataset is pretty accessible by the standards of modern datasets;
the code as-is runs an epoch in about ten seconds on my laptop's CPU (four on a Colab GPU).
Here we only run three epochs out of laziness, but if you train longer it'll do better.

\begin{qlist}
\item \pts{15}
Implement your own graph convolution operation, rather than using PyTorch Geometric's; there's scaffolding for you in \texttt{MyGCNConv}. (\texttt{main.py my-gcn} runs that for you; you shouldn't need to change anything in the \texttt{GCN} class.) A graph convolution is given by
\[
    \operatorname{GCN}(x)_v
    = \sum_{u \in \operatorname{neighbours}(v) \cup \{v\}} \frac{1}{\sqrt{\operatorname{deg}(v) \operatorname{deg}(u)}} W x_u
    + b
,\]
where $v$ and $u$ are nodes in the graph with corresponding node feature vectors $x_v$ and $x_u$.
(Recall that a neighbour of a node is any node with an edge  the degree of a node is its number of neighbours.)
The parameters of the matrix are the matrix $W$ of shape \texttt{[dim\_out, dim\_in]}
and the bias vector $b$ of shape \texttt{[dim\_out]}.

While they implement it with a relatively complex message passing framework, do not use this in your code.
You may want to use a matrix multiplication framing (you probably don't want to do explicit looping).
Your implementation will likely be slower than theirs, but run it for a few epochs to make sure the performance is about the same;
my straightforward matrix-multiplication implementation takes about 90 seconds per epoch (instead of 10) on my laptop's CPU,
or about 7 seconds (instead of 4) on a Colab GPU.
The accuracy might not be exactly the same (there's a bunch of randomness, and a few things will be slightly different);
as long as it's going down and not a huge amount higher than the other implementation, it should be fine.
\ask{Hand in your code.}

\begin{answer}\TODO\end{answer}

\item \pts{25}
While there are many variants of Transformers for graph data,
the class \texttt{gnns.GraphTransformer} implements a simple variant that does pairwise attention over all nodes in the graph.
Similarly to how the sequence order of a sentence is only available to a standard Transformer via positional encodings,
in this version of a graph Transformer, the structure of the graph is only available via the Laplacian positional encoding.
While this can cause issues on some datasets, in this dataset it turns out to be pretty okay.
\texttt{main.py graph-transformer} runs this model using torch's implementation of multi-headed self-attention.

Because of the additional pairwise attention,
this method is slower than GCNs,
so if you don't have a decent GPU it'll really be nicer to train on Colab (where an epoch takes 15 seconds, instead of about 10 minutes on my CPU!).
But you can develop locally, make sure your code runs for a couple of batches, and then try training on Colab.
(If you want, try training for longer to get a much better predictive model!)

\ask{Finish the implementation of the \texttt{gnns.MultiHeadSelfAttention.forward()} method};
you can run it with \texttt{main.py my-graph-transformer}.
My implementation is again slower than the pytorch implementation (which has had a lot of optimization effort put into it),
but regression performance is the same.
Make sure it runs okay, and \ask{hand in your code}.

\begin{answer}\TODO\end{answer}

\end{qlist}

\end{document}

