Bug fix in backward likelihood calculation (#142)

* use stationary(obj.model) in backward likelihood calculation, handle llg0 and llg1 = 0.0 in optimizegamma

* added comment block to avoid calculating oldlik and warnings after NNIs

Co-authored-by: Cecile Ane <cecileane@users.noreply.github.com>

src/phyLiNCoptimization.jl

@@ -1703,19 +1703,17 @@ function optimizegamma_LiNC!(obj::SSM, focusedge::Edge,
     # derivative of loglik at γ=0
     llg0 = wsum(ulik)
     inside01 = true
-    if llg0 < 0
+    if llg0 <= 0.0 # if llg0 = 0, something fishy is happening or data has no variation
         γ = 0.0
         inside01 = false
         ll = wsum((visib ? log.(clikp) : log.(clikp .+ clikn)))
         @debug "γ = 0 best, skip Newton-Raphson"
-    elseif llg0 == 0.0 # if llg0 = 0, something fishy is happening.
-        @error("llg0 is $llg0 and optimization is proceeding. Something is wrong.")
     else # at γ=1
         if visib
              ulik .= (clike .- clikp) ./ clike
         else ulik .= (clike .- clikp) ./ (clike .+ clikn); end
         llg1 = wsum(ulik)
-        if llg1 > 0
+        if llg1 >= 0.0 # if llg1 = 0.0, data has no variation
             γ = 1.0
             inside01 = false
             ll = wsum((visib ? log.(clike) : log.(clike .+ clikn)))

src/traitsLikDiscrete.jl

@@ -1074,7 +1074,7 @@ function discrete_backwardlikelihood_trait!(obj::SSM, t::Integer, ri::Integer)
     k = nstates(obj.model)
     fill!(backwardlik, 0.0) # re-initialize for each trait, each iteration
     bkwtmp = Vector{Float64}(undef, k) # to hold bkw lik without parent edge transition
-    if typeof(obj.model) == NASM
+    if typeof(obj.model) <: NASM
         logprior = log.(stationary(obj.model))
     else #trait models
         logprior = [-log(k) for i in 1:k] # uniform prior at root