\documentclass[notitlepage]{article}
\usepackage{../Math323}
\usepackage{listings}
\usepackage{numprint}
\usepackage{enumerate}
\usepackage{multirow}
\usepackage{bbm}

\usepackage{amsfonts,amsmath}

\setlength{\parindent}{0pt}

\def\E{\Expect}


\lstloadlanguages{R}

\definecolor{keywordcolor}{rgb}{0,0.6,0.6}
\definecolor{delimcolor}{rgb}{0.461,0.039,0.102}
\definecolor{Rcommentcolor}{rgb}{0.101,0.043,0.432}

\lstdefinestyle{Rsettings}{
  basicstyle=\ttfamily,
  breaklines=true,
  showstringspaces=false,
  keywords={if, else, function, theFunction, tmp}, % Write as many keywords
  otherkeywords={},
  commentstyle=\itshape\color{Rcommentcolor},
  keywordstyle=\color{keywordcolor},
  moredelim=[s][\color{delimcolor}]{"}{"},
}

\lstset{basicstyle=\ttfamily, numbers=none, literate={~} {$\sim$}{2}}

\begin{document}

<<setup, cache=FALSE, echo=FALSE>>=
library(knitr)
# global chunk options
opts_chunk$set(cache=TRUE, autodep=TRUE)
options(scipen=999)
options(repos=c(CRAN="https://cloud.r-project.org/"))
inline_hook <- function (x) {
  if (is.numeric(x)) {
    # ifelse does a vectorized comparison
    # If integer, print without decimal; otherwise print two places
    res <- ifelse(x == round(x),
      sprintf("%.3f", x),
      sprintf("%.3f", x)
    )
    paste(res, collapse = ", ")
  }
}
knit_hooks$set(inline = inline_hook)
@
\begin{center}
{\textsclarge{Joint Discrete Distributions}}
\end{center}

\textbf{The joint pmf:} The \textit{joint probability mass function} (joint pmf) for two discrete random variables $Y_1$ and $Y_2$ is the function
\[
p(y_1,y_2) \equiv p_{Y_1,Y_2}(y_1,y_2) = P(Y_1 = y_1, Y_2 = y_2)  \qquad (y_1,y_2) \in \Z \times \Z.
\]
We can tabulate $p_{Y_1,Y_2}(y_1,y_2)$ as follows:  suppose $Y_1$ takes values on $\{0,1,2,\ldots,8\}$ and $Y_2$ takes values on $\{0,1,2,\ldots,6\}$, then we may write
\begin{table}[ht]
\centering
\begin{tabular}{ll|ccccccccc|}
\hline
& & & & & & & & & &\\[-6pt]
\multirow{9}{*}{$Y_2$} & 6 &  $p(0,6)$ & $p(1,6)$ & $p(2,6)$ & $p(3,6)$ & $p(4,6)$ & $p(5,6)$& $p(6,6)$ & $p(7,6)$ & $p(8,6)$\\[3pt]
& 5 &  $p(0,5)$ & $p(1,5)$ & $p(2,5)$ & $p(3,5)$ & $p(4,5)$ & $p(5,5)$& $p(6,5)$ & $p(7,5)$ & $p(8,5)$\\[3pt]
& 4 &  $p(0,4)$ & $p(1,4)$ & $p(2,4)$ & $p(3,4)$ & $p(4,4)$ & $p(5,4)$& $p(6,4)$ & $p(7,4)$ & $p(8,4)$\\[3pt]
& 3 &  $p(0,3)$ & $p(1,3)$ & $p(2,3)$ & $p(3,3)$ & $p(4,3)$ & $p(5,3)$& $p(6,3)$ & $p(7,3)$ & $p(8,3)$\\[3pt]
& 2 &  $p(0,2)$ & $p(1,2)$ & $p(2,2)$ & $p(3,2)$ & $p(4,2)$ & $p(5,2)$& $p(6,2)$ & $p(7,2)$ & $p(8,2)$\\[3pt]
& 1 &  $p(0,1)$ & $p(1,1)$ & $p(2,1)$ & $p(3,1)$ & $p(4,1)$ & $p(5,1)$& $p(6,1)$ & $p(7,1)$ & $p(8,1)$\\[3pt]
& 0 &  $p(0,0)$ & $p(1,0)$ & $p(2,0)$ & $p(3,0)$ & $p(4,0)$ & $p(5,0)$& $p(6,0)$ & $p(7,0)$ & $p(8,0)$\\[3pt]
& & & & & & & & &\\[-12pt]
\hline
& & & & & & & & & &\\[-6pt]
& & 0 & 1 & 2 & 3 & 4 & 5 & 6 & 7 & 8\\[3pt]
& & \multicolumn{9}{c|}{$Y_1$}\\
\hline
\end{tabular}
\end{table}
For example, suppose that
\[
p(y_1,y_2) = c \exp\{1+(y_1-4)(y_2-2)/10\} \qquad 0 \leq y_1 \leq 8, 0 \leq y_2 \leq 6
\]
and $p(y_1,y_2) = 0$ otherwise, for some constant $c > 0$.  The constant $c$ ensures that the joint pmf sums to 1, that is
\[
c \sum_{y_1=0}^8 \sum_{y_2 = 0}^6 \exp\{1+(y_1-4)(y_2-2)/10\} = 1
\]
<<dists01,comment='+', fig.width=7.25, fig.height=4.00, fig.align='center', echo=FALSE>>=
pfunc<-function(x,y){return(exp((1 + (x-4)*(y-2))/10))}
ijmat<-t(outer(0:8,0:6,'+'))
ijmat<-ijmat[nrow(ijmat):1,]
#ijmat
p<-t(outer(0:8,0:6,pfunc))
cval<-1/sum(p)
pmat<-cval*p
pmat<-pmat[nrow(pmat):1,]
pmat<-t(pmat)
#round(pmat,3)
@
In this case we can compute $c=\Sexpr{round(cval,4)}$, and compute the probabilities as follows:
\begin{table}[ht]
\centering
\begin{tabular}{ll|ccccccccc|}
\hline
& & & & & & & & & &\\[-6pt]
\multirow{9}{*}{$Y_2$} & 6 &  $\Sexpr{round(pfunc(0,6)*cval,3)}$ & $\Sexpr{round(pfunc(1,6)*cval,3)}$ & $\Sexpr{round(pfunc(2,6)*cval,3)}$ & $\Sexpr{round(pfunc(3,6)*cval,3)}$ & $\Sexpr{round(pfunc(4,6)*cval,3)}$ & $\Sexpr{round(pfunc(5,6)*cval,3)}$& $\Sexpr{round(pfunc(6,6)*cval,3)}$ & $\Sexpr{round(pfunc(7,6)*cval,3)}$ & $\Sexpr{round(pfunc(8,6)*cval,3)}$\\[3pt]
& 5 &  $\Sexpr{round(pfunc(0,5)*cval,3)}$ & $\Sexpr{round(pfunc(1,5)*cval,3)}$ & $\Sexpr{round(pfunc(2,5)*cval,3)}$ & $\Sexpr{round(pfunc(3,5)*cval,3)}$ & $\Sexpr{round(pfunc(4,5)*cval,3)}$ & $\Sexpr{round(pfunc(5,5)*cval,3)}$& $\Sexpr{round(pfunc(6,5)*cval,3)}$ & $\Sexpr{round(pfunc(7,5)*cval,3)}$ & $\Sexpr{round(pfunc(8,5)*cval,3)}$\\[3pt]
& 4 &  $\Sexpr{round(pfunc(0,4)*cval,3)}$ & $\Sexpr{round(pfunc(1,4)*cval,3)}$ & $\Sexpr{round(pfunc(2,4)*cval,3)}$ & $\Sexpr{round(pfunc(3,4)*cval,3)}$ & $\Sexpr{round(pfunc(4,4)*cval,3)}$ & $\Sexpr{round(pfunc(5,4)*cval,3)}$& $\Sexpr{round(pfunc(6,4)*cval,3)}$ & $\Sexpr{round(pfunc(7,4)*cval,3)}$ & $\Sexpr{round(pfunc(8,4)*cval,3)}$\\[3pt]
& 3 &  $\Sexpr{round(pfunc(0,3)*cval,3)}$ & $\Sexpr{round(pfunc(1,3)*cval,3)}$ & $\Sexpr{round(pfunc(2,3)*cval,3)}$ & $\Sexpr{round(pfunc(3,3)*cval,3)}$ & $\Sexpr{round(pfunc(4,3)*cval,3)}$ & $\Sexpr{round(pfunc(5,3)*cval,3)}$& $\Sexpr{round(pfunc(6,3)*cval,3)}$ & $\Sexpr{round(pfunc(7,3)*cval,3)}$ & $\Sexpr{round(pfunc(8,3)*cval,3)}$\\[3pt]
& 2 &  $\Sexpr{round(pfunc(0,2)*cval,3)}$ & $\Sexpr{round(pfunc(1,2)*cval,3)}$ & $\Sexpr{round(pfunc(2,2)*cval,3)}$ & $\Sexpr{round(pfunc(3,2)*cval,3)}$ & $\Sexpr{round(pfunc(4,2)*cval,3)}$ & $\Sexpr{round(pfunc(5,2)*cval,3)}$& $\Sexpr{round(pfunc(6,2)*cval,3)}$ & $\Sexpr{round(pfunc(7,2)*cval,3)}$ & $\Sexpr{round(pfunc(8,2)*cval,3)}$\\[3pt]
& 1 &  $\Sexpr{round(pfunc(0,1)*cval,3)}$ & $\Sexpr{round(pfunc(1,1)*cval,3)}$ & $\Sexpr{round(pfunc(2,1)*cval,3)}$ & $\Sexpr{round(pfunc(3,1)*cval,3)}$ & $\Sexpr{round(pfunc(4,1)*cval,3)}$ & $\Sexpr{round(pfunc(5,1)*cval,3)}$& $\Sexpr{round(pfunc(6,1)*cval,3)}$ & $\Sexpr{round(pfunc(7,1)*cval,3)}$ & $\Sexpr{round(pfunc(8,1)*cval,3)}$\\[3pt]
& 0 &  $\Sexpr{round(pfunc(0,0)*cval,3)}$ & $\Sexpr{round(pfunc(1,0)*cval,3)}$ & $\Sexpr{round(pfunc(2,0)*cval,3)}$ & $\Sexpr{round(pfunc(3,0)*cval,3)}$ & $\Sexpr{round(pfunc(4,0)*cval,3)}$ & $\Sexpr{round(pfunc(5,0)*cval,3)}$& $\Sexpr{round(pfunc(6,0)*cval,3)}$ & $\Sexpr{round(pfunc(7,0)*cval,3)}$ & $\Sexpr{round(pfunc(8,0)*cval,3)}$\\[3pt]
& & & & & & & & &\\[-12pt]
\hline
& & & & & & & & & &\\[-6pt]
& & 0 & 1 & 2 & 3 & 4 & 5 & 6 & 7 & 8\\[3pt]
& & \multicolumn{9}{c|}{$Y_1$}\\
\hline
\end{tabular}
\end{table}

Thus for example
\[
p(0,3)=\Sexpr{round(pfunc(0,3)*cval,3)} \qquad p(4,1)=\Sexpr{round(pfunc(4,1)*cval,3)} \qquad p(3,6)=\Sexpr{round(pfunc(3,6)*cval,3)} \qquad p(8,5)=\Sexpr{round(pfunc(8,5)*cval,3)}
\]

This function specifies how much probability mass is placed at each $(y_1,y_2)$ value; the total probability mass is 1.  The function can be displayed in several different ways:

\begin{itemize}

\item \textbf{Surface plot:}  here the probability masses are plotted vertically at the specified $(y_1,y_2)$ positions, with a vertical line projecting down onto the $(y_1,y_2)$-plane.


<<dists02,comment='+', fig.width=8, fig.height=6, fig.align='center', echo=FALSE>>=
library(scatterplot3d)
x<-y<-z<-numeric()
for(y1 in 0:8){
    for(y2 in 0:6){
        x<-c(x,y1)
        y<-c(y,y2)
        z<-c(z,pfunc(y1,y2)*cval)
     }
}
jpmf<-data.frame(x,y,z)
scatterplot3d(jpmf[,1:3], pch = 16, type="h",color='red',ylim=range(-0.5,6.5),xlim=range(-0.5,8.5),zlim=range(0,0.1),scale.y=2/3,
xlab=expression(y[1]),ylab=expression(y[2]),zlab=expression(p(y[1],y[2])))
@

\item \textbf{Point plot:}  here the radius of the point at $(y_1,y_2)$ indicates the amount of probability placed there.

<<dists03,comment='+', fig.width=16/3, fig.height=4, fig.align='center', echo=FALSE>>=
par(mar=c(4,4,1,1))
plot(0:8,0:8,ylim=range(-0.5,6.5),xlim=range(-0.5,8.5),xlab=expression(y[1]),ylab=expression(y[2]),type='n',xaxt="n")
axis(1, at = seq(0, 8, by = 1))
for(y1 in 0:8){
    for(y2 in 0:6){
        points(y1,y2,cex=pfunc(y1,y2)*cval/max(pmat)*3,pch=19)
    }
}
@

\pagebreak

\item \textbf{Image plot:} here the colours represent the probability value (yellow is high, blue is low).

<<dists04,comment='+', fig.width=16/3, fig.height=4, fig.align='center', echo=FALSE>>=
par(mar=c(4,4,1,1))
xv<-seq(-0.5,8.5,by=1)
yv<-seq(-0.5,6.5,by=1)
pmat<-pmat[,ncol(pmat):1]
colfunc <- colorRampPalette(c("blue","red","yellow"))
image(xv,yv,pmat,xlab=expression(y[1]),ylab=expression(y[2]),xaxt="n",col=colfunc(100))
axis(1, at = seq(0, 8, by = 1))
@
\end{itemize}

From the joint pmf, we can compute the \textit{marginal} pmf for each of the variables $Y_1$ and $Y_2$.  We obtain these using the Theorem of Total Probability as follows: we have for example that
\begin{align*}
p_{Y_1}(y_1) = P(Y_1 = y_1) & = P(Y_1 = y_1, Y_2 = 0) + P(Y_1 = y_1, Y_2 = 1) + \cdots + P(Y_1 = y_1, Y_2 = 6)\\[6pt]
& = p_{Y_1,Y_2}(y_1,0) + p_{Y_1,Y_2}(y_1,1) + \cdots + p_{Y_1,Y_2}(y_1,6)
\end{align*}
and so
\begin{align*}
Y_1: \quad p_{Y_1}(y_1) & = \sum_{y_2=0}^6 p_{Y_1,Y_2}(y_1,y_2)\\[6pt]
Y_2: \quad p_{Y_2}(y_2) & = \sum_{y_1=0}^8 p_{Y_1,Y_2}(y_1,y_2).
\end{align*}
Thus we obtain the marginal distribution for $Y_1$ by summing down the \textit{columns} of the above tables, and the marginal distribution for $Y_2$ by summing across the \textit{rows} of the tables.  In the above numerical example, we have
<<dists05,comment='+', fig.width=8, fig.height=6, fig.align='center', echo=FALSE>>=
p1<-rep(0,9)
p2<-rep(0,7)
for(y1 in 0:8){
    for(y2 in 0:6){
       p1[y1+1]<-p1[y1+1]+pfunc(y1,y2)*cval
       p2[y2+1]<-p2[y2+1]+pfunc(y1,y2)*cval
    }
}
@
\begin{table}[ht]
\centering
\begin{tabular}{|l|ccccccccc|}
\hline
$y_1$ & 0 & 1 & 2 & 3 & 4 & 5 & 6 & 7 & 8\\[3pt]
$p_{Y_1}(y_1)$ &  $\Sexpr{round(p1[1],3)}$ & $\Sexpr{round(p1[2],3)}$ & $\Sexpr{round(p1[3],3)}$ & $\Sexpr{round(p1[4],3)}$ & $\Sexpr{round(p1[5],3)}$ & $\Sexpr{round(p1[6],3)}$ & $\Sexpr{round(p1[7],3)}$ & $\Sexpr{round(p1[8],3)}$ & $\Sexpr{round(p1[9],3)}$\\
\hline
\end{tabular}
\end{table}
and
\begin{table}[ht]
\centering
\begin{tabular}{|l|ccccccc|}
\hline
$y_2$ & 0 & 1 & 2 & 3 & 4 & 5 & 6 \\[3pt]
$p_{Y_2}(y_2)$ &  $\Sexpr{round(p2[1],3)}$ & $\Sexpr{round(p2[2],3)}$ & $\Sexpr{round(p2[3],3)}$ & $\Sexpr{round(p2[4],3)}$ & $\Sexpr{round(p2[5],3)}$ & $\Sexpr{round(p2[6],3)}$ & $\Sexpr{round(p2[7],3)}$ \\
\hline
\end{tabular}
\end{table}

Note that these functions are both pmfs: we have
\[
\sum_{y_1=0}^8 p_{Y_1}(y_1) = \sum_{y_2=0}^6 p_{Y_2}(y_2) = 1.
\]

\pagebreak
We may plot these pmfs in the usual way.

<<dists06,comment='+', fig.width=8, fig.height=4, fig.align='center', echo=FALSE>>=
par(mar=c(4,4,4,1),mfrow=c(1,2))
plot(0:8,p1,pch=19,xlab=expression(y[1]),ylab=expression(p[Y[1]](y[1])),ylim=range(0,0.3))
for(i in 0:8){lines(c(i,i),c(0,p1[i+1]))}
title(expression(p[Y[1]](y[1])))
axis(1, at = seq(0, 8, by = 1))
plot(0:6,p2,pch=19,xlab=expression(y[2]),ylab=expression(p[Y[2]](y[2])),ylim=range(0,0.3))
for(i in 0:6){lines(c(i,i),c(0,p2[i+1]))}
title(expression(p[Y[2]](y[2])))
axis(1, at = seq(0, 8, by = 1))
@

Finally, we can consider \textit{conditional} pmfs:  we have by the definition of conditional probability
\[
P(Y_1 = y_1 | Y_2 = y_2) = \frac{P(Y_1 = y_1 , Y_2 = y_2)}{P(Y_2 = y_2)} \qquad 0 \leq y_1 \leq 8
\]
for each $y_2$, $0 \leq y_2 \leq 6$.  This enables us to define the conditional pmf for $Y_1$, given that $Y_2 = y_2$, as
\[
p_{Y_1|Y_2}(y_1|y_2) = P(Y_1 = y_1 | Y_2 = y_2) = \frac{p_{Y_1,Y_2}(y_1,y_2)}{p_{Y_2}(y_2)} \qquad 0 \leq y_1 \leq 8
\]
for each $y_2$, $0 \leq y_2 \leq 6$.  This function is a pmf in terms of $y_1$, that is for example
\[
\sum_{y_1=0}^8 p_{Y_1|Y_2}(y_1|y_2) = 1
\]
for each $y_2$, $0 \leq y_2 \leq 6$.     Similarly we may define the conditional pmf for $Y_2$, given that $Y_1 = y_1$, as
\[
p_{Y_2|Y_1}(y_2|y_1) = P(Y_2 = y_2 | Y_1 = y_1) = \frac{p_{Y_1,Y_2}(y_1,y_2)}{p_{Y_1}(y_1)} \qquad 0 \leq y_2 \leq 6
\]
for each $y_1$, $0 \leq y_1 \leq 8$.  This function is a pmf in terms of $y_2$.

\medskip
We compute the conditional pmfs by taking `slices' through the joint pmf, that is, we examine individual rows or columns, and then standardize by dividing by the row or column sum.   For example, for the conditional pmf for $Y_1$ given that $Y_2=4$, we examine the $Y_2=4$ row
\begin{table}[ht]
\centering
\begin{tabular}{|ll|ccccccccc|}
\hline
& & & & & & & & & &\\[-6pt]
$Y_4$ & 4 &  $p(0,4)$ & $p(1,4)$ & $p(2,4)$ & $p(3,4)$ & $p(4,4)$ & $p(5,4)$& $p(6,4)$ & $p(7,4)$ & $p(8,4)$\\[3pt]
& & & & & & & & &\\[-12pt]
\hline
& & & & & & & & & &\\[-6pt]
& & 0 & 1 & 2 & 3 & 4 & 5 & 6 & 7 & 8\\[3pt]
& & \multicolumn{9}{c|}{$Y_1$}\\
\hline
\end{tabular}
\end{table}

and divide by $p_{Y_2}(4)$ to get the conditional pmf.

\pagebreak

Numerically, we have for $p_{Y_1|Y_2}(y_1|4)$
\begin{table}[ht]
\centering
\begin{tabular}{|ll|ccccccccc|}
\hline
& & & & & & & & & &\\[-6pt]
$Y_2$ & 4 &  \Sexpr{round(pfunc(0,4)*cval/p2[5],3)} & \Sexpr{round(pfunc(1,4)*cval/p2[5],3)} & \Sexpr{round(pfunc(2,4)*cval/p2[5],3)} & \Sexpr{round(pfunc(3,4)*cval/p2[5],3)} & \Sexpr{round(pfunc(4,4)*cval/p2[5],3)} & \Sexpr{round(pfunc(5,4)*cval/p2[5],3)} & \Sexpr{round(pfunc(6,4)*cval/p2[5],3)} & \Sexpr{round(pfunc(7,4)*cval/p2[5],3)} & \Sexpr{round(pfunc(8,4)*cval/p2[5],3)}\\[3pt]
& & & & & & & & &\\[-12pt]
\hline
& & & & & & & & & &\\[-6pt]
& & 0 & 1 & 2 & 3 & 4 & 5 & 6 & 7 & 8\\[3pt]
& & \multicolumn{9}{c|}{$Y_1$}\\
\hline
\end{tabular}
\end{table}


<<dists07,comment='+', fig.width=8, fig.height=5, fig.align='center',echo=FALSE>>=
par(mar=c(4,4,4,1))
p1.4<-p1.2<-rep(0,9)
for(y1 in 0:8){
    p1.2[y1+1]<-pfunc(y1,2)*cval
    p1.4[y1+1]<-pfunc(y1,4)*cval
}
p1.2<-p1.2/sum(p1.2)
p1.4<-p1.4/sum(p1.4)
plot((0:8)-0.1,p1.2,pch=19,xlab=expression(y[1]),ylab='Prob',xlim=range(-0.5,8.5),ylim=range(0,0.25),col='blue')
for(i in 0:8){lines(c(i-0.1,i-0.1),c(0,p1.2[i+1]),col='blue')}
axis(1, at = seq(0, 8, by = 1))
points((0:8)+0.1,p1.4,pch=19,col='red')
for(i in 0:8){lines(c(i+0.1,i+0.1),c(0,p1.4[i+1]),col='red')}
title(expression(p[Y[1]*"|"*Y[2]](y[1]*"|"*y[2])))
legend(0,0.25,c(expression(p[Y[1]*"|"*Y[2]](y[1]*"|"*4)),expression(p[Y[1]*"|"*Y[2]](y[1]*"|"*2))),col=c('red','blue'),pch=19)
@

For the conditional pmf for $Y_1$ given that $Y_2=2$, we examine the $Y_2=2$ row
\begin{table}[ht]
\centering
\begin{tabular}{|ll|ccccccccc|}
\hline
& & & & & & & & & &\\[-6pt]
$Y_2$ & 2 &  $p(0,2)$ & $p(1,2)$ & $p(2,2)$ & $p(3,2)$ & $p(4,2)$ & $p(5,2)$& $p(6,2)$ & $p(7,2)$ & $p(8,2)$\\[3pt]
& & & & & & & & &\\[-12pt]
\hline
& & & & & & & & & &\\[-6pt]
& & 0 & 1 & 2 & 3 & 4 & 5 & 6 & 7 & 8\\[3pt]
& & \multicolumn{9}{c|}{$Y_1$}\\
\hline
\end{tabular}
\end{table}
and divide by $p_{Y_2}(2)$ to get the conditional pmf.    Numerically, we have for $p_{Y_1|Y_2}(y_1|2)$
\begin{table}[ht]
\centering
\begin{tabular}{|ll|ccccccccc|}
\hline
& & & & & & & & & &\\[-6pt]
$Y_2$ & 2 &  \Sexpr{round(pfunc(0,2)*cval/p2[3],3)} & \Sexpr{round(pfunc(1,2)*cval/p2[3],3)} & \Sexpr{round(pfunc(2,2)*cval/p2[3],3)} & \Sexpr{round(pfunc(3,2)*cval/p2[3],3)} & \Sexpr{round(pfunc(4,2)*cval/p2[3],3)} & \Sexpr{round(pfunc(5,2)*cval/p2[3],3)} & \Sexpr{round(pfunc(6,2)*cval/p2[3],3)} & \Sexpr{round(pfunc(7,2)*cval/p2[3],3)} & \Sexpr{round(pfunc(8,2)*cval/p2[3],3)}\\[3pt]
& & & & & & & & &\\[-12pt]
\hline
& & & & & & & & & &\\[-6pt]
& & 0 & 1 & 2 & 3 & 4 & 5 & 6 & 7 & 8\\[3pt]
& & \multicolumn{9}{c|}{$Y_1$}\\
\hline
\end{tabular}
\end{table}

that is, constant across all values of $y_1$, $0 \leq y_1 \leq 8$.  In general, the function $p_{Y_1|Y_2}(y_1|y_2)$ is different for each $y_2$,

\pagebreak

Similarly, to get the conditional pmf for $Y_2$ given that $Y_1=5$, we examine the $Y_1=5$ column and divide by $p_{Y_1}(5)$ to get the conditional pmf
\begin{table}[ht]
\centering
\begin{tabular}{|ll|ccccccc|}
\hline
& & & & &  & & &\\[-6pt]
$Y_1$ & 5 &  \Sexpr{round(pfunc(5,0)*cval/p2[3],3)} & \Sexpr{round(pfunc(5,1)*cval/p2[3],3)} & \Sexpr{round(pfunc(5,2)*cval/p2[3],3)} & \Sexpr{round(pfunc(5,3)*cval/p2[3],3)} & \Sexpr{round(pfunc(5,4)*cval/p2[3],3)} & \Sexpr{round(pfunc(5,5)*cval/p2[3],3)} & \Sexpr{round(pfunc(5,6)*cval/p2[3],3)} \\[3pt]
& & & & & & &\\[-12pt]
\hline
& & & & & & & &\\[-6pt]
& & 0 & 1 & 2 & 3 & 4 & 5 & 6 \\[3pt]
& & \multicolumn{7}{c|}{$Y_2$}\\
\hline
\end{tabular}
\end{table}

To get the conditional pmf for $Y_2$ given that $Y_1=0$, we examine the $Y_1=0$ column and divide by $p_{Y_1}(0)$ to get the conditional pmf
\begin{table}[ht]
\centering
\begin{tabular}{|ll|ccccccc|}
\hline
& & & & &  & & &\\[-6pt]
$Y_1$ & 5 &  \Sexpr{round(pfunc(0,0)*cval/p1[1],3)} & \Sexpr{round(pfunc(0,1)*cval/p1[1],3)} & \Sexpr{round(pfunc(0,2)*cval/p1[1],3)} & \Sexpr{round(pfunc(0,3)*cval/p1[1],3)} & \Sexpr{round(pfunc(0,4)*cval/p1[1],3)} & \Sexpr{round(pfunc(0,5)*cval/p1[1],3)} & \Sexpr{round(pfunc(0,6)*cval/p1[1],3)} \\[3pt]
& & & & & & &\\[-12pt]
\hline
& & & & & & & &\\[-6pt]
& & 0 & 1 & 2 & 3 & 4 & 5 & 6 \\[3pt]
& & \multicolumn{7}{c|}{$Y_2$}\\
\hline
\end{tabular}
\end{table}

<<dists08,comment='+', fig.width=8, fig.height=5, fig.align='center',echo=FALSE>>=
par(mar=c(4,4,4,1))
p2.0<-p2.5<-rep(0,7)
for(y2 in 0:6){
    p2.0[y2+1]<-pfunc(0,y2)*cval
    p2.5[y2+1]<-pfunc(5,y2)*cval
}
p2.0<-p2.0/sum(p2.0)
p2.5<-p2.5/sum(p2.5)
plot((0:6)-0.1,p2.0,pch=19,xlab=expression(y[2]),ylab='Prob',xlim=range(-0.5,6.5),ylim=range(0,0.5),col='blue')
for(i in 0:6){lines(c(i-0.1,i-0.1),c(0,p2.0[i+1]),col='blue')}
axis(1, at = seq(0, 6, by = 1))
points((0:6)+0.1,p2.5,pch=19,col='red')
for(i in 0:6){lines(c(i+0.1,i+0.1),c(0,p2.5[i+1]),col='red')}
title(expression(p[Y[2]*"|"*Y[1]](y[2]*"|"*y[1])))
legend(0,0.5,c(expression(p[Y[2]*"|"*Y[1]](y[2]*"|"*0)),expression(p[Y[2]*"|"*Y[1]](y[2]*"|"*5))),col=c('red','blue'),pch=19)
@

\end{document}