Hi - I have a data set (~900 x 120) on which I’m trying to do dimension reduction. Eventually, I want to test some factor mixture models, but right now, I simply want to do a basic factor analysis. The data set isn’t *huge*, but it’s large enough that I’m not expecting things to be fast. However, it’s even slower than I expect.

Here’s a MWE to give a sense of things, on 1.4.2, Windows 64-bit:

```
using Turing, Optim
data1 = randn(900, 120)
# p = the dimension of the factor model
# this is just an example; we ignore rotational invariance issues
@model factormodel(dat, p) = begin
N, K = size(dat)
σ ~ filldist(InverseGamma(1, 1), K) # factor uniquenesses
λ ~ filldist(Normal(0, 1), p, K) # factor loadings
μ ~ filldist(Normal(0, 10), K) # variable-specific intercept
F ~ filldist(Normal(0,1), N, p) # latent factors
mu = F * λ
for k in 1:K
for i in 1:N
dat[i,k] ~ Normal(μ[k] + mu[i,k], σ[k])
end
end
end;
#=
(1) MLE/MAP estimates take an *extraordinarily* long time, especially compared to e.g., R's factanal() or fa() routines. I would expect it to take a few seconds, maybe. Instead, it takes so long I get tired of waiting and cancel. :( I'm not willing to wait overnight for a model like this, ha ha.
(2) standard NUTS estimation also takes a *very* long time (and seems to spend a lot of time generating values that are "rejected for numerical errors").
(3) basic MH is slightly slow, but usable (-ish).
(4) Going from ForwardDiff to a reverse-mode AD (like Zygote) does not change anything.
=#
# Trying a reverse diff backend like Zygote doesn't change anything.
# using Zygote
# Turing.setadbackend(:zygote)
map_estimate = optimize(factormodel(data1, 1), MAP())
hmc_estimate = sample(factormodel(data1, 1), NUTS(0.65), 2000)
mh_estimate = sample(factormodel(data1, 1), MH(), 20000)
# Also tested ADVI on master, but gives errors
# advi_estimate = vi(factormodel(data1, 1), ADVI(10, 1000))
```

I tested on both the current stable Turing (0.13.0) and master, with same results. I certainly don’t expect this kind of model to be *fast*, but I would expect it to run in a feasible time. A simple one-factor model only has 900+240 parameters, which is a good number, but not outlandish.

Is there something obvious I’m missing or doing very wrong? I have read the manual’s section on optimization, and I know the one loop isn’t ideal, but I’m not sure what to replace it with, or if it’s the real culprit here. I’d welcome any advice people have.