Allow constructing log-scaled distributions directly

[bredelings]

This allows constructing distributions like logNormal by supplying the distribution in log space:

x ~ dnLog(dnNormal(0,1))
x ~ dnLog(dnGamma(a,b))

Currently the only log-scaled distributions implemented are dnLognormal, dnLoguniform, and dnLogExponential. This PR allows us to construct (for example)

dnLog(dnGamma(a,b)), dnLog(dnLaplace(m,s)), dnLog(dnCauchy(m,s)), etc.

And equally importantly, future distributions that have not been implemented yet.

[hoehna]

Overall this looks very good and sensible. I have only one theoretical question: Is there a Jaccobian needed? In other words, do you get the same between those cases:

x1 ~ dnExp(5)
y1 := ln( x1 )
moves.append( mvScale(x1) )

y2 ~ dnLogExponential(5)
moves.append( mvSlide(y2) )

y3 ~ dnLog( dnExp(5) )
moves.append( mvSlide(y3) )

All three should be equivalent, right?

[bredelings]

Yes, a Jacobian is included. If we have dnLog(f), then the density is f(log(x))/x. The 1/x term is the Jacobian.

However, my understanding is that x ~ logExponential(5) means log(x) ~ Exponential(5). What you wrote here is the opposite though:

x1 ~ dnExp(5)
y1 := ln( x1 )
moves.append( mvScale(x1) )

Since x ~ logNormal(mu,sigma)) means that log(x) ~ Normal(mu,sigma), I think we need to do it the other way. It is unfortunate that a logNormal sample is not the log of a Normal sample.

I checked that the density for dnLog(dnNormal(0,1)) that matches the density for dnLognormal(0,1):

> y1 ~ dnLognormal(0,1)
> y1.clamp(2)
> y1.probability()
   0.156874
> y2 ~ dnLog(dnNormal(0,1))
> y2.clamp(2)
> y2.probability()
   0.156874
> y2.clamp(3)
> y2.probability()

The density for dnLogUniform(exp(-1),exp(1)) doesn't match the density for dnLog(dnUniform(-1,1)) or the theoretical density. But it looks like it is only off by a constant factor of log(10), because it uses log10 instead of log to compute the density.

I notice that RevBayes does not match R in that RevBayes uses ln(x) for the natural log, whereas R allows a base but the base defaults to exp(1).

[bredelings]

Hmm... I notice that the help files for dnLoguniform and dnLogExponential both claim that they are uniform on orders of magnitude. It cannot be true for both of them. I can fix up dnLogExponential so that it has the correct pdf, but then its not clear what its purpose is. So perhaps we should delete it and people can use dnLogUniform instead. Any code that was using dnLogExponential was already not working correctly, I think.

[hoehna]

Why did you remove the distribution? This breaks older scripts that used this distribution. One option could be to internally create LogDistribution if you want to remove some c++ files. The same should apply to the LogUniform distribution?

[bredelings]

My reasoning is that there is no option here that is without problems. I think dnLogExponential had the wrong density before. At first I fixed the density, but that silently changes the meaning of the older scripts, while allowing them to run. For example, previously x ~ dnLogExponential(1) could return a negative x. However, when the correct density is used, x is always positive. So it seems like the previous scripts need to be modified somehow. This option is not great because it breaks the older scripts -- but at least it tells the user how to fix them.

I don't like breaking the older scripts, but I also worry about changing their meaning silently. I see that dnLogExponential is used only in one tutorial that was written by Mike May in 2022. We could ask him what he wants to do...

[bredelings]

OK, now you can also write:

Mu = [1.0, 2.0, 3.0, 4.0]
Sigma ~ dnWishart(3,2,4)
x ~ dnLog(dnMultivariateNormal(Mu,Sigma))