initial_values seem to fail with LKJ priors

Hi! Whenever the `lkj_correlation()` prior is part of a model, the `initial_values` argument seems to always cause `mcmc()` to fail.

Reproducible code from the `greta` [example models page](https://greta-stats.org/articles/example_models.html#ecological-models):
```
# model matrix
modmat <- model.matrix(~ Sepal.Width, iris)
# index of species
jj <- as.numeric(iris$Species)

M <- ncol(modmat) # number of varying coefficients
N <- max(jj) # number of species

# prior on the standard deviation of the varying coefficient
tau <- exponential(0.5, dim = M)

# prior on the correlation between the varying coefficient
Omega <- lkj_correlation(3, M)

# optimization of the varying coefficient sampling through
# cholesky factorization and whitening
Omega_U <- chol(Omega)
Sigma_U <- sweep(Omega_U, 2, tau, "*")
z <- normal(0, 1, dim = c(N, M))
ab <- z %*% Sigma_U # equivalent to: ab ~ multi_normal(0, Sigma_U)

# the linear predictor
mu <- rowSums(ab[jj,] * modmat)

# the residual variance
sigma_e <- cauchy(0, 3, truncation = c(0, Inf))

#model
y <- iris$Sepal.Length
distribution(y) <- normal(mu, sigma_e)
m <- model(ab, sigma_e)

draws <- mcmc(m, chains = 4, initial_values = initials(sigma_e = 1))
```

On my computer this throws an error:
```
Error in py_call_impl(callable, dots$args, dots$keywords) :
ValueError: Cannot feed value of shape (1, 13) for Tensor 'Placeholder:0', which has shape '(?, 10)'

Detailed traceback:
File "/home/hrlai/.local/share/r-miniconda/envs/r-reticulate/lib/python3.6/site-packages/tensorflow/python/client/session.py", line 950, in run
run_metadata_ptr)
File "/home/hrlai/.local/share/r-miniconda/envs/r-reticulate/lib/python3.6/site-packages/tensorflow/python/client/session.py", line 1149, in _run
str(subfeed_t.get_shape())))
```

On various datasets, the `Cannot feed value of shape (1, 13) for Tensor 'Placeholder:0', which has shape '(?, 10)'` message always have the first number `13` larger than the second number `10`. And the difference between them seems to be `N` --- this led me to suspect the LKJ prior in the first place.

I tried to remove LKJ from the model and initial values work again. (After restarting R session) This is reproducible via:
```
# model matrix
modmat <- model.matrix(~ Sepal.Width, iris)
# index of species
jj <- as.numeric(iris$Species)

M <- ncol(modmat) # number of varying coefficients
N <- max(jj) # number of species

# prior on the standard deviation of the varying coefficient
tau <- exponential(0.5, dim = M)

Sigma_U <- zeros(dim = c(M, M))
diag(Sigma_U) <- tau
z <- normal(0, 1, dim = c(N, M))
ab <- z %*% Sigma_U # equivalent to: ab ~ multi_normal(0, Sigma_U)

# the linear predictor
mu <- rowSums(ab[jj,] * modmat)

# the residual variance
sigma_e <- cauchy(0, 3, truncation = c(0, Inf))

#model
y <- iris$Sepal.Length
distribution(y) <- normal(mu, sigma_e)
m <- model(ab, sigma_e)

draws <- mcmc(m, chains = 4, initial_values = initials(sigma_e = 1))
```

I'd really like to keep the LKJ prior as well as being able to specify initial values to help chain convergence. Looking forward to hear your idea!