Add beta, remove beta location scale and add te (#120)

Author: guilhermebodin

Project.toml

@@ -1,7 +1,7 @@
 name = "ScoreDrivenModels"
 uuid = "4a87933e-d659-11e9-0e65-7f40dedd4a3a"
 authors = ["guilhermebodin <[email protected]>, raphaelsaavedra <[email protected]>"]
-version = "0.1.4"
+version = "0.1.5"
 
 [deps]
 Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f"

docs/src/manual.md

@@ -95,7 +95,6 @@ this is handled internally.
 | Distribution | Identity scaling | Inverse and inverse square root scalings |
 | :------- |:---------------:|:--------------------:|
 |[`Beta`](@ref)| ✓ | ✓  |
-|[`BetaLocationScale`](@ref)| ✓ |  x  |
 |[`Exponential`](@ref)| ✓ | ✓  |
 |[`Gamma`](@ref)| ✓ | ✓  |
 |[`LogitNormal`](@ref)| ✓ | ✓  |
@@ -110,7 +109,6 @@ this is handled internally.
 
 ```@docs
 ScoreDrivenModels.Beta
-ScoreDrivenModels.BetaLocationScale
 ScoreDrivenModels.Exponential
 ScoreDrivenModels.Gamma
 ScoreDrivenModels.LogitNormal

src/MLE.jl

@@ -94,7 +94,8 @@ end
 function fit(gas::Model{D, T}, y::Vector{T};
              initial_params::Matrix{T} = DEFAULT_INITIAL_PARAM,
              opt_method::AbstractOptimizationMethod = NelderMead(gas, DEFAULT_NUM_SEEDS),
-             verbose::Int = DEFAULT_VERBOSE) where {D, T}
+             verbose::Int = DEFAULT_VERBOSE,
+             throw_errors::Bool = false) where {D, T}
 
     verbose in [0, 1, 2, 3] || throw(ErrorException, "verbose argument must be in [0, 1, 2, 3]")
     # Number of initial_points and number of params to estimate
@@ -123,6 +124,9 @@ function fit(gas::Model{D, T}, y::Vector{T};
             verbose >= 1 && println("Round $i of $n_initial_points - Log-likelihood: $(-opt_result.minimum * length(y))")
         catch err
             println(err)
+            if throw_errors
+                throw(err)
+            end
             verbose >= 1 && println("Round $i diverged")
         end
     end

src/ScoreDrivenModels.jl

@@ -31,7 +31,7 @@ include("opt_methods/NelderMead.jl")
 # Distributions
 include("distributions/non_native_dists.jl")
 include("distributions/common_interface.jl")
-include("distributions/betalocationscale.jl")
+include("distributions/beta.jl")
 include("distributions/exponential.jl")
 include("distributions/gamma.jl")
 include("distributions/logitnormal.jl")

src/distributions/beta.jl

@@ -0,0 +1,63 @@
+"""
+    Beta
+* Parametrization
+* Score
+* Fisher Information
+* `time_varying_params` map.
+* Default link
+"""
+Beta
+
+function score!(score_til::Matrix{T}, y::T, ::Type{Beta}, param::Matrix{T}, t::Int) where T
+    score_til[t, 1] = log(y) + digamma(param[t, 1] + param[t, 2]) - digamma(param[t, 1])
+    score_til[t, 2] = log(1 - y) + digamma(param[t, 1] + param[t, 2]) - digamma(param[t, 2])
+    return
+end
+
+function fisher_information!(aux::AuxiliaryLinAlg{T}, ::Type{Beta}, param::Matrix{T}, t::Int) where T
+    minus_trigamma_a_b = -trigamma(param[t, 1] + param[t, 2])
+    aux.fisher[1, 1] = trigamma(param[t, 1]) + minus_trigamma_a_b
+    aux.fisher[2, 2] = trigamma(param[t, 2]) + minus_trigamma_a_b
+    aux.fisher[2, 1] = minus_trigamma_a_b
+    aux.fisher[1, 2] = minus_trigamma_a_b
+    return
+end
+
+function log_likelihood(::Type{Beta}, y::Vector{T}, param::Matrix{T}, n::Int) where T
+    loglik = 0.0
+    for t in 1:n
+        loglik += (param[t, 1] - 1)*log(y[t]) + (param[t, 2] - 1)*log(1 - y[t]) - logbeta(param[t, 1], param[t, 2])
+    end
+    return -loglik
+end
+
+# Links
+function link!(param_tilde::Matrix{T}, ::Type{Beta}, param::Matrix{T}, t::Int) where T 
+    param_tilde[t, 1] = link(LogLink, param[t, 1], zero(T))
+    param_tilde[t, 2] = link(LogLink, param[t, 2], zero(T))
+    return
+end
+function unlink!(param::Matrix{T}, ::Type{Beta}, param_tilde::Matrix{T}, t::Int) where T 
+    param[t, 1] = unlink(LogLink, param_tilde[t, 1], zero(T))
+    param[t, 2] = unlink(LogLink, param_tilde[t, 2], zero(T))
+    return
+end
+function jacobian_link!(aux::AuxiliaryLinAlg{T}, ::Type{Beta}, param::Matrix{T}, t::Int) where T 
+    aux.jac[1] = jacobian_link(LogLink, param[t, 1], zero(T))
+    aux.jac[2] = jacobian_link(LogLink, param[t, 2], zero(T))
+    return
+end
+
+# utils
+function update_dist(::Type{Beta}, param::Matrix{T}, t::Int) where T
+    small_threshold!(param, SMALL_NUM, t)
+    return Beta(param[t, 1], param[t, 2])
+end 
+
+function params_sdm(d::Beta)
+    return Distributions.params(d)
+end
+
+function num_params(::Type{Beta})
+    return 2
+end

src/distributions/betalocationscale.jl

@@ -1,6 +1,6 @@
 # Currently supported distributions
 const DISTS = [
-    BetaLocationScale,
+    Beta,
     Exponential,
     Gamma,
     LogitNormal,
@@ -13,7 +13,7 @@ const DISTS = [
     Weibull
 ]
 
-export BetaLocationScale,
+export Beta,
     Exponential,
     Gamma,
     LogitNormal,

src/distributions/common_interface.jl

@@ -1,7 +1,3 @@
 # Define TDistLocationScale as the location scale transformation of the Distributions.jl TDist
 export TDistLocationScale
-TDistLocationScale = LocationScale{Float64,TDist{Float64}}
-
-# Define BetaLocationScale as the location scale transformation of the Distributions.jl TDist
-export BetaLocationScale
-BetaLocationScale = LocationScale{Float64,Beta{Float64}}
+TDistLocationScale = LocationScale{Float64,TDist{Float64}}

src/distributions/non_native_dists.jl

@@ -1,5 +1,5 @@
 const PARAMETER_REQUIRED_DISTS = [Chisq; Chi; TDist]
-const FI_NOT_IMPLEMENTED = [Weibull; BetaLocationScale; NegativeBinomial]
+const FI_NOT_IMPLEMENTED = [Weibull; NegativeBinomial]
 
 struct FakeDist{T<:Real} <: Distributions.ContinuousUnivariateDistribution
     foo::T
@@ -9,9 +9,6 @@ function instantiate_dist(D::Type{<:Distribution})
     if D == TDistLocationScale
         params = [1.0 2.0 10]
         return ScoreDrivenModels.update_dist(D, params, 1)
-    elseif D == BetaLocationScale
-        params = [-2.0 5.0 1.5 3.4]
-        return ScoreDrivenModels.update_dist(D, params, 1)
     elseif D in PARAMETER_REQUIRED_DISTS
         params = [10][:, :]
         return ScoreDrivenModels.update_dist(D, params, 1)