This gives GPflow likelihoods a stronger contract re what input shapes are expected and what shapes are returned. In particular, we should obey something akin to tensorflow_probability’s event-shape/batch-shape/sample-shape. Very little changes for most users, except that some shapes will be asserted. Advanced users will benefit from more shape checks and better defined return shapes for methods attached to likelihoods.

Likelihoods now need to define:
 - self.observation_dim: what is the last dimension of Y supposed to be?
 - self.latent_dim: what is the last dimension of F, F_mu and F_var expected to be?
We’ll check that the dimensions of tensors passed in match up.
Return shapes for all methods will be the broadcast-shape of the tensors passed in, with the last dimension removed.
Example: likelihood.variational_expectations(F_mu, F_var, Y) might take tensors of dimensions [..., 2], [..., 2], [..., 2] and return a tensor of shape [...]

The shape checks are handled by the public methods log_prob, predict_mean_and_var, predict_log_density, and variational_expectations; new likelihoods should implement the corresponding private methods with leading underscore.

## Standard likelihoods
Most likelihoods in GPflow are univariate, and treat columns of Y and F as independent variables. For these observation_dim and latent_dim are None, and we shape-check F and Y on-the-fly for matching last-column dimensions. Note that the return shape contract changes as per that above, and the likelihood methods return the sum over observation dimensions.

## Fancy likelihoods
Likelihoods that depart from the univariate standard include:
SwitchedLikelihood [we’ll check that latent_dim = observation_dim - 1]
MultiClass/Softmax [observation_dim = 1, latent_dim = number of classes (e.g. 10 for MNIST)]
HeteroskedasticGaussian e.g. see GPflow notebook [observation_dim = 1, latent_dim = 2]

Note that this change deprecates Likelihood.predict_density in favour of Likelihood.predict_log_density.