Introduction to Neural Networks
U. Minn. Psy 5038

Gaussian generative models, learning, and inference

Initialize standard library files:

Off[General::spell1];

<<"BarCharts`";

Last time

Quick review of probability and statistics
Applications to random sampling

Today

Review of big picture
Examples of computations on continuous probabilities
Examples of computations on discrete probabilities
Introduction to Bayes learning

Recall relationship between neural networks and Bayesian inference

Here we are only talking about inference or estimation based on patterns of neural activity--i.e. in the language of neural networks, about "recall", rather than learning. Later we will introduce Bayesian learning. 

In the  general case we can talk about the probability over all possible values of a neural network's state vector: p((Subscript[V, 1]), (Subscript[V, 2]),...)

This doesn't distinguish which values are fixed and which are allowed to vary. If some values are fixed, then we can treat those as the input, and allow the network's free neurons to vary to maximize a conditional probability:

Relationship between posterior probability and the energy of the state of a Hopfield neural network:

H (the hypothesis space) corresponds to the values of state variables, i.e. patterns of neural activity that are changing {(Subscript[V, 1]), (Subscript[V, 2]),...}. d corresponds to "data" or the fixed, clamped values {(Subscript[V, 1]^s), (Subscript[V, 2]^s)...}

From the more abstract point of view of statistical inference, we have some variables that are fixed--the "data" d, some we let vary to maximize probability--the hypotheses H, and some that we may not care to estimate (n). Let's see how to use these distinctions.
Suppose we have p(H,d,n). Two rules of inference:

Marginalize over what you don't care about:

p(H,d)=(∑_n)p(H,d,n)

Condition on what you know:

p(H|d)=((p(H,d))/(p(d)))

Mathematica functions for gaussian multivariates & exploring marginals

Define PDF, CDF

Needs["MultivariateStatistics`"]

The add-on function, ({{MultinormalDistribution[μ,Σ]  specifies a, "multinormal (multivariate Gaussian) distribution with mean vector "μ" and covariance matrix "Σ" "}})

m1 = {1,1/2};
r=(1/2)*{{1,2/3},{2/3,4}};
ndist = MultinormalDistribution[m1, r];
pdf = PDF[ndist, {x1, x2}]

g1=ContourPlot[PDF[ndist,{x1,x2}],{x1,-3,3},{x2,-3,3}]

What is the probability of the distribution in the region (Subscript[x, 1])<.5⋂(Subscript[x, 2])<2?

grp=RegionPlot[x1<.5 && x2<2,{x1,-4,4},{x2,-4,4},PlotStyle->Directive[Opacity[.25],EdgeForm[],FaceForm[Gray]]];

Show[{ g1,grp},ImageSize->Small]

gcdf = ContourPlot[CDF[ndist, {x1, x2}], {x1, -4, 4}, {x2, -4, 4},ImageSize -> Small];
Show[{ gcdf, grp}, ImageSize -> Small]

CDF[ndist, {.5, 2.0}]

Clear[x1,x2];
marginal[x1_]:=(∫_(-∞)^∞)PDF[ndist,{x1,x2}]\[DifferentialD]x2;
marginal2[x2_]:=(∫_(-∞)^∞)PDF[ndist,{x1,x2}]\[DifferentialD]x1;

mt=Table[{x1,marginal[x1]},{x1,-3,3,.2}];

g2=ListPlot[mt,Joined->True,PlotStyle->{Red,Thick},Axes->False];

mt2=Table[{x2,marginal2[x2]},{x2,-3,3,.4}];
g3=ListPlot[mt2,Joined->True,PlotStyle->{Green,Thick},Axes->False];

theta=Pi/2;
Show[g1,Epilog->{Inset[g2,{0,-3},{0,0}],Inset[g3,{-3,0},{0,0},Automatic,{{Cos[theta],Sin[theta]},{Sin[theta],-Cos[theta]}}]}]

Finding the mode

FindMaximum[PDF[ndist,{x1,x2}],{{x1,0},{x2,0}}]

Drawing samples

As we've used in earlier lectures, drawing samples is done by:

RandomReal[ndist]

Mixtures of gaussians  with MultinormalDistribution[]

Clear[mix];

r1=0.4*{{1,.6},{.6,1}};
r2=0.4*{{1,-.6},{-.6,1}};
m1 = {1,.5}; m2 = {-1,-.5};
ndist1 = MultinormalDistribution[m1, r1];
ndist2 = MultinormalDistribution[m2, r2];

mix[x_] := 0.5 (PDF[ndist1, x] + PDF[ndist2, x]);

gg1=ContourPlot[mix[{x1,x2}],{x1,-2,2},{x2,-2,2},PlotRange->Full,ImageSize->Small]

Marginals for mixture

marginal[x1_]:=Integrate[mix[{x1,x2}],{x2,-Infinity,Infinity}]

Clear[marginal];
marginal[x1_]:=0.5*(NIntegrate[PDF[ndist1,{x1,x2}],{x2,-Infinity,Infinity}] + NIntegrate[PDF[ndist2,{x1,x2}],{x2,-Infinity,Infinity}]);

gg2=Plot[marginal[x1],{x1,-2,2},PlotStyle->{Red,Thick},Axes->{False,False}];

Clear[marginal];
marginal[x2_]:=0.5*(NIntegrate[PDF[ndist1,{x1,x2}],{x1,-Infinity,Infinity}] + NIntegrate[PDF[ndist2,{x1,x2}],{x1,-Infinity,Infinity}]);

gg3=Plot[0.0333*marginal[x2],{x2,-2,2},PlotStyle->{Green,Thick},Axes->{False,False}];

theta=Pi/2;
Show[gg1,ImageSize->200,Epilog->{Inset[gg2,{0,-3},{0,0}],Inset[gg3,{-2,0},{0,0},Automatic,{{Cos[theta],Sin[theta]},{Sin[theta],-Cos[theta]}}]}]

Which projection (marginal) is more "interesting"--the one onto x1 or onto x2? 
Exploratory projection pursuit

Graphical Models of dependence

Graphs: causal structure and conditional independence

In general, natural image pattern formation is specified by a high-dimensional joint probability, requiring an elaboration of the causal structure that is more complex than the simple SDT model.  The idea is to represent the probabilistic structure of the joint distribution P(S,L,I) by a Bayes net (e.g. Ripley, 1996}, which is a graphical model that expresses how variables influence each other.  There are  three basic building blocks: converging, diverging, and intermediate nodes.  For example, multiple (e.g. scene) variables causing a given image measurement, a single variable producing multiple image measurements, or a cause indirectly influencing an image measurement through an intermediate variable.  These types of influence provide a first step towards modeling the joint distribution and the means to compute probabilities of the unknown variables given known values.

Components of the generative structure for image patterns involve converging, diverging,and intermediate nodes. For example,these could correspond to:multiple (scene) causes {shape S1, illumination S2 giving rise to the same image measurement, I ; one cause, S influencing more than one image measurement, {color, I1, brightness, I2}; a scene (or other) cause S, {object identity, S} influencing an image measurement (image contour) through an intermediate variable L (3D shape) .

The arrows tell us how to factor the joint probability into conditionals. So for the three examples above, we have:
p(S1,S2,I)=p(I|S1,S2)p(S1)p(S2)
p(S,I1,I2)=p(I1|S)p(I2|S)p(S)
p(S,L,I)=p(I|L)p(L|S)p(S)

Primary, secondary variables.

The following figure draws a parallel between the causal structure, as determined by the generative model, for signal detection theory (as in the light detection problem), and the general problem of visual inference.

We can interpret the causal structure in terms of conditional probability.
The top panel shows one possible generative graph structure for an  ideal observer problem in classical signal detection theory (SDT).  The data are determined by the signal hypotheses plus (additive gaussian) noise.  Knowledge is represented by the joint probability p(x,H,n)=p(x|H,n)p(H)p(n).  The lower panel shows a simplified example of the generative structure for perceptual inference from a pattern inference theory perspective.  The image measurements (x) are determined by a typically non-linear function (\phi) of primary signal variables (S_e) and confounding secondary variables (S_g).  Knowledge is represented by the joint probability p(x,S_e,S_g).  Both scene and image variables can be high  dimensional vectors. In general, the causal structure of natural image patterns is more complex and consequently requires elaboration of its graphical representation. For SDT and pattern inference theory, the task is to make a decision about the signal hypotheses or primary signal variables, while discounting the noise or secondary variables.  Thus optimal perceptual decisions are determined by p(x,S_e), which is derived by summing over the secondary variables (i.e. marginalizing with respect to the secondary variables): (∫_(S_g))p(x,S_e,S_g)dS_g.

Influences between variables are represented by conditioning, and a graphical model expresses the conditional independencies between variables.  Two random variables may only become independent, however, once the value of some third variable is known.  This is called conditional independence.Recall from above that two random variables are independent if and only if their joint probability is equal to the product of their individual probabilities.  Thus, if p(A,B) = p(A)p(B), then A and B are independent.  If p(A,B|C) = p(A|C)p(B|C), then A and B are conditionally independent.  When corn prices drop in the summer, hay fever incidence goes up.  However, if the joint on corn price and hay fever is conditioned on ``ideal weather for corn and ragweed'', the correlation between corn prices and hay fever drops.  This is because Corn price and hay fever symptoms are conditionally independent. 
There is a correlation between eating ice cream and drowning. Why? What event should you condition on to make the dependence go away?

What is noise? Primary and secondary variables in SDT and in pattern inference theory

Noise is whatever you don't care to estimate, but contributes to the data.

Optimal Inference and task dependence: Fruit example

(due to James Coughlan; see Yuille, Coughlan, Kersten & Schrater).

Figure from Yuille, Coughlan, Kersten & Schrater.

The graph specifies how to decompose the joint probability:
 p[F, C, Is, Ic ] = p[ Ic | C ] p[C | F ] p[Is | F ] p[F ]

Generative model: The prior model on hypotheses, F & C

More apples (F=1) than tomatoes (F=2), and:

ppF[F_]:=If[F==1,9/16,7/16];
TableForm[Table[ppF[F],{F,1,2}],TableHeadings->{{"F=a", "F=t"}}]

The conditional  probability cpCF[C|F]:

cpCF[F_,C_]:=Which[F==1 && C ==1,5/9,F==1 && C ==2,4/9,F==2 && C ==1,6/7,F==2 && C ==2,1/7];
TableForm[Table[cpCF[F,C],{C,1,2},{F,1,2}],TableHeadings->{{"C=r","C=g"},{"F=a", "F=t"}}]

So the joint is:

jpFC[F_,C_]:=cpCF[F,C] ppF[F];
TableForm[Table[jpFC[F,C],{F,1,2},{C,1,2}],TableHeadings->{{"F=a", "F=t"},{"C=r","C=g"}}]

We can marginalize to get the prior probability on color alone is:

ppC[C_]:=(∑_(F=1)^2)jpFC[F,C]

Question: Is fruit identity independent of material color--i.e. is F independent of C? Check whether the joint probability on Fruit and Color can be factored into the product of the prior probabilities on Fruit and Color.

Answer

No.

TableForm[Table[jpFC[F,C],{F,1,2},{C,1,2}],TableHeadings->{{"F=a", "F=t"},{"C=r","C=g"}}]
TableForm[Table[ppF[F] ppC[C],{F,1,2},{C,1,2}],TableHeadings->{{"F=a", "F=t"},{"C=r","C=g"}}]

Generative model: The likelihood model--probabilities of measurements, i.e. some features given hypotheses

Suppose that we have gathered some "image statistics" which provides us knowledge of how the image measurements for shape Is, and for color Ic depend on the type of fruit F, and material color, C. For simplicity, our measurements are discrete and binary (a more realistic case, they would have continuous values), say Is = {am, tm}, and Ic = {rm, gm}.

P(I_S= am, tm | F=a) = {11/16, 5/16}
P(I_S= am, tm | F=t) = {5/8, 3/8}
P(I_C= rm, gm | C=r) = {9/16, 7/16}
P(I_C= rm, gm | C=g) = {1/2, 1/2}
We use the notation am, tm, rm, gm because the measurements are already suggestive of the likely cause. So there is a correlation between apple and apple-like shapes, am; and between red material, and "red" measurements. In general, there may not be an obvious correlation like this.

We define a function for the  probability of Ic given C,  cpIcC[Ic | C]:

cpIcC[Ic_,C_]:=Which[Ic==1 && C ==1,9/16,Ic==1 && C ==2,7/16,Ic==2 && C ==1,1/2,Ic==2 && C ==2,1/2];
TableForm[Table[cpIcC[Ic,C],{C,1,2},{Ic,1,2}],TableHeadings->{{"Ic=rm","Ic=gm"},{"C=r", "C=g"}}]

The  probability of Is conditional on F is cpIsF[Is | F]:

cpIsF[Is_,F_]:=Which[Is==1 && F ==1,11/16,Is==1 && F ==2,5/8,Is==2 && F ==1,5/16,Is==2 && F ==2,3/8];
TableForm[Table[cpIsF[Is,F],{Is,1,2},{F,1,2}],TableHeadings->{{"Is=am","Is=tm"},{"F=a", "F=t"}}]

The total joint probability

We now have enough information to put probabilities on the 2x2x2 "universe" of possibilities, i.e. all possible combinations of fruit, color, and image measurements. Looking at the graphical model makes it easy to use the product rule to construct the total joint, which is:
 p[F, C, Is, Ic ] = p[ Ic | C ] p[C | F ] p[Is | F ] p[F ]:

jpFCIsIc[F_, C_, Is_, Ic_ ] := cpIcC[ Ic, C ] cpCF[F,C ] cpIsF[Is,F ] ppF[F ]

Usually, we don't need the probabilities of the image measurements (because once the measurements are made, they are fixed and we want to compare the probabilities of the hypotheses. But in our simple case here, once we have the joint, we can calculate the probabilities of the image measurements through marginalization p(Is,Ic)=(∑_C)(∑_F)p(F,C,Is,Ic), too:

jpIsIc[Is_,Ic_]:=(∑_(C=1)^2)(∑_(F=1)^2)jpFCIsIc[F, C, Is, Ic ]

Three MAP tasks

We are going to show that the best guess (i.e. maximum probability) depends on the task.

Define argmax[] function:

argmax[x_]:=Position[x,Max[x]];

Pick most probable fruit AND color--Answer "red tomato"

First, suppose the task is to make the best bet as to the fruit AND material color. To make it concrete, suppose that we see an "apple-ish shape" with a reddish color, i.e., we measure Is=am=1, and Ic = rm=1. The measurements suggest "red apple", but to find the most probable, we need to take into account the priors too in order to make the best guesses. p(F,C | Is, Ic) is given by:

FCcIsIc[F_,C_,Is_,Ic_]:=jpFCIsIc[F, C, Is, Ic ]/jpIsIc[Is,Ic]

TableForm[FCcIsIcTable=Table[FCcIsIc[F, C, 1, 1 ],{F,1,2},{C,1,2}],TableHeadings->{{"F=a", "F=t"},{"C=r","C=g"}}]
Max[FCcIsIcTable]
argmax[FCcIsIcTable]

Note that we don't have to calculate the conditional. We can use the fact that p(F,C | Is, Ic) = ((p(F,C,Is,Ic))/(p(Is,Ic)))∝ p(F, C , Is=1, Ic=1). I.e. the conditional is proportional to the function specified by the joint evaluated at Is=1, and Ic=1.

TableForm[jpFCIsIcTable=Table[jpFCIsIc[F, C, 1, 1 ],{F,1,2},{C,1,2}],TableHeadings->{{"F=a", "F=t"},{"C=r","C=g"}}]
Max[jpFCIsIcTable]
argmax[jpFCIsIcTable]

Question: Is p(F, C, Is=1, Ic=1) a joint probability on F and C? (Answer: No.)

Total[jpFCcIsIcTable,2]
Total[jpFCIsIcTable,2]

In either case, we conclude that "Red tomato" is the most probable once we take into account the difference in priors.

Pick most probable color--Answer "red"

Same measurements as before. But now suppose we only care about the true material color, and not the identity of the object. Then we want to integrate out or marginalize with respect to the shape or fruit-type variable, F. In this case, we want to maximize the posterior:

p(C | Is=1, Ic=1)=(∑_(F=1)^2)p(F,C|Is=1,Ic=1)

pC[C_,Is_,Ic_]:=(∑_(F=1)^2)jpFCcIsIc[F, C, Is, Ic ]

TableForm[pCTable=Table[pC[C, 1, 1 ],{C,1,2}],TableHeadings->{{"C=r","C=g"}}]
Max[pCTable]
argmax[pCTable]

Answer is that the most probable material color is C = r, "red".

Pick most probable fruit--Answer "apple"

Same measurements as before. But now, we don't care about the material color, just the identity of the fruit. Then we want to integrate out or marginalize with respect to the material variable, C. In this case, we want to maximize the posterior:

p(F | Is, Ic)

pF[F_,Is_,Ic_]:=(∑_(C=1)^2)jpFCcIsIc[F, C, Is, Ic ]

TableForm[pFTable=Table[pF[F, 1, 1 ],{F,1,2}],TableHeadings->{{"F=a","F=t"}}]
Max[pFTable]
argmax[pFTable]

The answer is "apple". So to sum up, for the same data measurements, the most probable fruit AND color is "red tomato", but the most probable fruit is "apple"!

Important "take-home message": Optimal inference depends on the precise definition of the task

Try expressing the consequences using the frequency interpretation of probability.

Bayesian learning of univariate Gaussian mean: MAP

From a statistical point of view, one form of learning is "density estimation" from histogram measurements. In high dimensions this is hard, but is easier if we have a low-dimensional parametric model for the density--i.e. the density is modeled in terms of a few parameters. So for example, the 1D Gaussian could be approximated by a  huge list of numbers--one for each bin, each number is an estimate of the probability of the value of the random variable falling in that bin. But because it is Gaussian, we can be more efficient by representing the density in terms of just two numbers (mean and variance), and a formula.
In this context, learning becomes parameter estimation.

A Bayesian learning example: Suppose we know the data comes from a Gaussian generative process, but we don't know the mean?

Suppose we have a set of samples that come from a Gaussian distribution with known variance (σ^2), but unknown mean μ.

(Subscript[x, i])=noise, where noise~N[μ,σ], or equivalently
(Subscript[x, i])=μ+noise, where noise~N[0,σ]

ndist0=NormalDistribution[μ,σ];

Although we don't know the mean, we can assume a Gaussian prior on the mean:

μ~N[μ0,σ0]

I.e. we make an initial guess of the mean's mean (μ0) and standard deviation (σ0). But we are willing to change our estimate of the mean given new data--i.e. given the posterior. If we are really uncertain at the beginning,, we can start of with a large standard deviation, and as we gather data, the uncertainty about the value of the mean will decrease.

ndistμ=NormalDistribution[μ0,(Subscript[σ, 0])];
PDF[ndistμ,μ]

Suppose the generative model N[μ,σ] produces three i.i.d. (independent, identically distributed) samples (Subscript[x, 1]),(Subscript[x, 2]),(Subscript[x, 3]).What is the MAP estimate of μ? Which value of  μ makes the posterior biggest? We use Bayes rule:

p(μ|(Subscript[x, 1]),(Subscript[x, 2]),(Subscript[x, 3]))=((p((Subscript[x, 1]),(Subscript[x, 2]),(Subscript[x, 3])|μ)p(μ))/(p((Subscript[x, 1]),(Subscript[x, 2]),(Subscript[x, 3]))))

p((Subscript[x, 1])|μ)is given by:

PDF[ndist0,(Subscript[x, 1])]

Because the samples are drawn independently, the p((Subscript[x, 1]),(Subscript[x, 2]),(Subscript[x, 3])|μ) is the product of three terms, so the numerator is p((Subscript[x, 1])|μ)p((Subscript[x, 2])|μ)p((Subscript[x, 3])|μ)times the prior p(μ):

PDF[ndist0,(Subscript[x, 1])]*PDF[ndist0,(Subscript[x, 2])]*PDF[ndist0,(Subscript[x, 3])]*PDF[ndistμ,μ]

Calculating the MAP estimate of mean

To find the value of the mean that is largest given our three samples, and our prior assumption, we find μ where p((Subscript[x, 1]),(Subscript[x, 2]),(Subscript[x, 3])|μ)p(μ)is biggest:

g=PDF[ndist0,(Subscript[x, 1])]*PDF[ndist0,(Subscript[x, 2])]*PDF[ndist0,(Subscript[x, 3])]*PDF[ndistμ,μ];
t=Log[g];
t1=PowerExpand[t]
t2=D[t1,μ]
Solve[-t2==0,μ]

In general, one can update from n samples in batch mode:

{{μ->((((μ0)/((σ0)^2))+(1/(σ^2))(∑_(i=1)^n)(Subscript[x, i]))/((n/(σ^2))+(1/((σ0)^2))))}}

For the multi-variate case, see Duda and Hart.

What is the influence of the initial estimate of the mean as learning goes on? What is the estimate of the mean as n gets large?

Appendices

Using Mathematica lists to manipulate discrete priors, likelihoods, and posteriors

A note on list arithmetic

We haven't done standard matrix/vector operations above to do conditioning. We've take advantage of how  Mathematica divides a 2x3 array by a 2-element vector:

M=Array[m,{2,3}]
X = Array[x,{2}]

M/X

Putting the probabilities back together again to get the joint

Transpose[Transpose[pHx] px]

pxH pH

Getting the posterior from the priors and likelihoods:

One reason Bayes' theorem is so useful is that it is often easier to formulate the likelihoods (e.g. from a causal or generativemodel of how the data could have occurred), and the priors (often from heuristics, or in computational vision empirically testable models of the external visual world). So let's use Mathematica to derive p(H|x) from p(x|H) and p(H) , (i.e. pHx from pxH and pH ).

px2 = Plus@@(pxH pH)

Transpose[Transpose[(pxH pH)]/Plus@@(pxH pH)]

Show that this joint probability has a uniform prior (i.e. both priors equal).

p={{1/8,1/8,1/4},{1/4,1/8,1/8}}

Marginalization and conditioning: A small dimensional example using list manipulation in Mathematica

A discrete joint probability

All of our knowledge regarding the signal discrimination problem can be described in terms of the joint probability of the hypotheses, H and the possible data measurements, x. The probability function assigns a number to all possible combinations:

p[H, x]

That is, we are assuming that both the hypotheses and the data are discrete random variables.

H={({{S1}, {S2}})

x ϵ{1,2,...}

Let's assume that x can only take on one of three values, 1, 2, or 3. And suppose the joint probability is:

p={{(1/12),(1/12),(1/6)},{(1/3),(1/6),(1/6)}}

TableForm[p,TableHeadings->{{"H=S1", "H=S2"},{"x=1","x=2","x=3"}}]

The total probability should sum up to one. Let's test to make sure. We first turn the list of lists into a singel list of scalars using Flatten[]. And then we can sum either with Apply[Plus,Flatten[p]].

Plus @@ Flatten[p]

We can pull out the first row of p like this:

p[[1]]

{(1/12),(1/12),(1/6)}

Is this the probability of x? No. For a start, the numbers don't sum to one. But we can get it through the two processes of marginalization and conditioning.

Marginalizing

What are the probabilities of the data, p(x)? To find out, we use the sum rule to sum over the columns:

px = Apply[Plus,p]

"Summing over "is also called marginalization or "integrating out".  Note that marginalization turns a probability function with higher degrees of freedom into one of lower degrees of freedom.

What are the prior probabilities? p(H)? To find out, we sum over the rows:

pH = Apply[Plus, Transpose[p]]

Conditioning

Now that we have the marginals, we can get use the product rule to obtain the conditional probability through conditioning of the joint:

p[ x|H]=((p[H, x])/(p[H]))

In the Exercises, you can see how to use Mathematica to do the division for conditioning. The syntax is simple:

pxH = p/pH

Note that the probability of x conditional on H sums up to 1 over x, i.e. each row adds up to 1. But, the columns do not. p[x|H] is a probability function of x, but a likelihood function of H. The posterior probability is obtained by conditioning on x:

p[ H|x]=((p[H, x])/(p[x]))

Syntax here is a bit more complicated, because the number of columns of px don't match the number of rows of p. We use Transpose[] to exchange the columns and rows of p before dividing, and then use Transpose again to get back the 2x3 form:

pHx = Transpose[Transpose[p]/px]

Plotting the joint

The following BarChart[] graphics function requires in add-in package (<< Graphics`Graphics`), which is specified at the top of the notebook. You could also use ListDensityPlot[].

BarChart[p[[1]],p[[2]]]

Marginalization and conditioning: An example using Mathematica functions

A discrete joint probability

All of our knowledge regarding the signal discrimination problem can be described in terms of the joint probability of the hypotheses, H and the possible data measurements, x. The probability function assigns a number to all possible combinations:

p[H, x]

That is, we are assuming that both the hypotheses and the data are discrete random variables.

H={({{S1}, {S2}})

x ϵ{1,2,...}

Let's assume that x can only take on one of three values, 1, 2, or 3. And suppose the joint probability is:

p[H_,x_]:=Which[H==1 && x ==1,1/12,H==1 && x ==2,1/12,H==1 && x ==3,1/6,H==2 && x ==1,1/3,H==2 && x ==2,1/6, H==2 && x ==3,1/6];

TableForm[Table[p[H,x],{H,1,2},{x,1,3}],TableHeadings->{{"H=s1", "H=s2"},{"X=1","X=2","X=3"}}]

The total probability should sum up to one. Let's test to make sure. We first turn the list of lists into a singel list of scalars using Flatten[]. And then we can sum either with Apply[Plus,Flatten[p]].

Sum[p[H,x],{H,1,2},{x,1,3}]

We can pull out the first row of p like this:

Table[p[1,x],{x,1,3}]

Is this the probability of x? No. For a start, the numbers don't sum to one. But we can get it through the two processes of marginalization and conditioning.

Marginalizing

What are the probabilities of the data, p(x)? To find out, we use the sum rule to sum over the columns:

px[x_]:=Sum[p[H,x],{H,1,2}];

Table[px[x],{x,1,3}]

"Summing over "is also called marginalization or "integrating out".  Note that marginalization turns a probability function with higher degrees of freedom into one of lower degrees of freedom.

What are the prior probabilities? p(H)? To find out, we sum over the rows:

pH[H_]:=Sum[p[H,x],{x,1,3}];

Table[pH[H],{H,1,2}]

Conditioning

Now that we have the marginals, we can get use the product rule to obtain the conditional probability through conditioning of the joint:

p[ x|H]=((p[H, x])/(p[H]))

We use function definition in Mathematica to do the division for conditioning. The syntax is simple:

pxH[H_,x_]:=p[H,x]/pH[H];

Table[pxH[H,x],{H,1,2},{x,1,3}]

Note that the probability of x conditional on H sums up to 1 over x, i.e. each row adds up to 1. But, the columns do not. p[x|H] is a probability function of x, but a likelihood function of H. The posterior probability is obtained by conditioning on x:

p[ H|x]=((p[H, x])/(p[x]))

pHx[H_,x_]:=p[H,x]/px[x];

Table[pHx [H,x],{H,1,2},{x,1,3}]

Plotting the joint

We use ListDensityPlot[].

ListDensityPlot[Table[p[H,x],{H,1,2},{x,1,3}]]

Random number generator, a non-Gaussian example: The von Mises distribution, with Matlab code

(courtesy, Paul Schrater)

function pofx = vonMisespdf(x,mu,sigma)
% For -pi <= x <= pi
% force x-mu within -pi to pi
y = angle(exp(i*(x-mu)));
kappa = 1/(sigma)^2;
%kappa = sigma;
pofx = exp(kappa*cos(y))/(2*pi*besseli(0,kappa));

function vonrand = vonMisesrand(nrand,mu,sigma)
% inverse cumulative method, executed by table lookup with
% linear interpolation
% build sampled cdf
x = (-pi:2*pi/(2e3):pi);
pofx = vonMisespdf(x,0,sigma);
cofx = cumsum(pofx/sum(pofx));
u = rand(1,nrand);
vonrand = interp1(cofx,x,u)+mu;

References

Applebaum, D. (1996). Probability and Information . Cambridge, UK: Cambridge University Press.
Cover, T. M., & Joy, A. T. (1991). Elements of Information Theory. New York: John Wiley & Sons, Inc.
Duda, R. O., & Hart, P. E. (1973). Pattern classification and scene  analysis . New York.: John Wiley & Sons.
Golden, R. (1988). A unified framework for connectionist systems. Biological Cybernetics, 59, 109-120.
Intrator, N.  Combining Exploratory Projection Pursuit and Projection Pursuit Regression.  Neural Computation (5):443-455, 1993. http://www.physics.brown.edu/users/faculty/intrator/papers/epp-ppr.ps.gz
Kersten, D. and P.W. Schrater (2000), Pattern Inference Theory: A Probabilistic Approach to Vision, in Perception and the Physical World, R. Mausfeld and D. Heyer, Editors. , John Wiley & Sons, Ltd.: Chichester. (pdf)
Kersten, D., Mamassian P & Yuille A (in press) Object perception as Bayesian inference. Annual Review of Psychology. (pdf, http://arjournals.annualreviews.org/doi/pdf/10.1146/annurev.psych.55.090902.142005)
Kersten, D., & Yuille, A. (2003). Bayesian models of object perception. Current Opinion in Neurobiology, 13(2) http://gandalf.psych.umn.edu/~kersten/kersten-lab/papers/KerstenYuilleApr2003.pdf
Ripley, B. D. (1996). Pattern Recognition and Neural Networks. Cambridge, UK: Cambridge University Press.
Yuille, A., Coughlan J., Kersten D.(1998) (pdf)