Merge branch 'main' into dev

NEWS.md

@@ -1,3 +1,13 @@
+# ModelConstructors.jl v0.2.4
+- Extend `transform_to_model_space` and `transform_to_real_line` for regime-switching parameters.
+
+# ModelConstructors.jl v0.2.3
+- Move methods for computing free and fixed indices of parameters from SMC.jl to ModelConstructors.jl
+
+# ModelConstructors.jl v0.2.2
+- Try to convert types to match rather than throwing a `MethodError` immediately
+  when calling `parameter`.
+
 # ModelConstructors.jl v0.2.1
 - Raise compat bounds for some packages
 

Project.toml

@@ -1,7 +1,7 @@
 name = "ModelConstructors"
 uuid = "e47e5152-bd14-11e9-1b46-c951f0a7041d"
 authors = ["William Chen <william.chen@ny.frb.org>", "Shlok Goyal <shlok.goyal@ny.frb.org>", "Alissa Johnson <alissa.johnson@ny.rb.org"]
-version = "0.2.1"
+version = "0.2.4"
 
 [deps]
 DataFrames = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0"

docs/src/index.md

@@ -30,5 +30,7 @@ Developers of this package at the
 [New York Fed](https://www.newyorkfed.org/research) include
 
 * [William Chen](https://github.com/chenwilliam77)
+* [Shlok Goyal](https://github.com/ShlokG)
+* [Alissa Johnson](https://github.com/alissarjohnson)
 * [Ethan Matlin](https://github.com/ethanmatlin)
 * [Reca Sarfati](https://github.com/rsarfati)

src/abstractmodel.jl

@@ -164,12 +164,66 @@ n_pseudo_observables(m::AbstractModel)      = length(m.pseudo_observables)
 n_equilibrium_conditions(m::AbstractModel)  = length(m.equilibrium_conditions)
 n_parameters(m::AbstractModel)              = length(m.parameters)
 n_parameters_steady_state(m::AbstractModel) = length(m.steady_state)
-n_parameters_free(m::AbstractModel)         = sum([!α.fixed for α in m.parameters])
+n_parameters_free(m::AbstractModel)         = length(get_free_para_inds(m.parameters; regime_switching = true))
 
 function n_parameters_regime_switching(m::AbstractModel)
     return n_parameters_regime_switching(m.parameters)
 end
 
+function get_fixed_para_inds(parameters::ParameterVector; regime_switching::Bool = false,
+                             toggle::Bool = true)
+    if regime_switching
+        if toggle
+            toggle_regime!(parameters, 1)
+        end
+
+        reg_fixed = [θ.fixed for θ in parameters] # it is assumed all regimes are toggled to regime 1
+        for θ in parameters
+            if !isempty(θ.regimes) # this parameter has regimes
+                if haskey(θ.regimes, :fixed)
+                    push!(reg_fixed, [regime_fixed(θ, i) for i in 2:length(θ.regimes[:value])]...)
+                elseif θ.fixed # since regimes[:fixed] is non-existent but θ.fixed is true,
+                    # it is assumed all regimes are fixed.
+                    push!(reg_fixed, trues(length(θ.regimes[:value]) - 1)...)
+                else # All regimes are not fixed
+                    push!(reg_fixed, falses(length(θ.regimes[:value]) - 1)...)
+                end
+            end
+        end
+
+        return findall(reg_fixed)
+    else
+        return findall([θ.fixed for θ in parameters])
+    end
+end
+
+function get_free_para_inds(parameters::ParameterVector; regime_switching::Bool = false,
+                            toggle::Bool = true)
+    if regime_switching
+        if toggle
+            toggle_regime!(parameters, 1)
+        end
+
+        reg_free = [!θ.fixed for θ in parameters] # it is assumed all regimes are toggled to regime 1
+        for θ in parameters
+            if !isempty(θ.regimes) # this parameter has regimes
+                if haskey(θ.regimes, :fixed)
+                    push!(reg_free, [!regime_fixed(θ, i) for i in 2:length(θ.regimes[:value])]...)
+                elseif θ.fixed # since regimes[:fixed] is non-existent but θ.fixed is true,
+                    # it is assumed all regimes are fixed.
+                    push!(reg_free, falses(length(θ.regimes[:value]) - 1)...)
+                else # All regimes are not fixed
+                    push!(reg_free, trues(length(θ.regimes[:value]) - 1)...)
+                end
+            end
+        end
+
+        return findall(reg_free)
+    else
+        return findall([!θ.fixed for θ in parameters])
+    end
+end
+
 """
 ```
 get_dict(m, class, index)

src/parameters.jl

@@ -387,7 +387,7 @@ function parameter(key::Symbol,
                    scaling::Function        = identity,
                    regimes::Dict{Symbol,OrderedDict{Int64,Any}} = Dict{Symbol,OrderedDict{Int64,Any}}(),
                    description::String = "No description available.",
-                   tex_label::String = "") where {V<:Vector, T <: Float64, U <:Transform} #{V<:Vector, S<:Real, T <: Float64, U <:Transform}
+                   tex_label::String = "") where {V<:Vector, T <: Real, U <:Transform} #{V<:Vector, S<:Real, T <: Float64, U <:Transform}
 
     # If fixed=true, force bounds to match and leave prior as null.  We need to define new
     # variable names here because of lexical scoping.
@@ -445,6 +445,32 @@ function parameter(key::Symbol,
     end
 end
 
+function parameter(key::Symbol,
+                   value::Union{T1, V}, #value::Union{S,V},
+                   valuebounds::Interval{T2} = (value,value),
+                   transform_parameterization::Interval{T3} = (value,value),
+                   transform::U             = Untransformed(),
+                   prior::Union{NullableOrPriorUnivariate, NullableOrPriorMultivariate} = NullablePriorUnivariate();
+                   fixed::Bool              = true,
+                   scaling::Function        = identity,
+                   regimes::Dict{Symbol,OrderedDict{Int64,Any}} = Dict{Symbol,OrderedDict{Int64,Any}}(),
+                   description::String = "No description available.",
+                   tex_label::String = "") where {V<:Vector, T1 <: Real, T2 <: Real, T3 <: Real, U <:Transform}
+    warn_str = "The element types of the fields `value` ($(typeof(value))), `valuebounds` ($(eltype(valuebounds))), " *
+        "and `transform_parameterization` ($(eltype(transform_parameterization))) do not match. " *
+        "Attempting to convert all types to the same type as `value`. Note that the element type for the prior " *
+        "distribution should also be $(typeof(value))."
+    @warn warn_str
+
+    valuebounds_new = (convert(T1, valuebounds[1]), convert(T1, valuebounds[2]))
+    transform_parameterization_new = (convert(T1, transform_parameterization[1]),
+                                      convert(T1, transform_parameterization[2]))
+
+    return parameter(key, value, valuebounds_new, transform_parameterization_new,
+                     transform, prior; fixed = fixed, scaling = scaling,
+                     regimes = regimes, description = description, tex_label = tex_label)
+end
+
 function parameter_ad(key::Symbol,
                       value::Union{S,V},
                       valuebounds::Interval{T} = (value,value),
@@ -750,8 +776,60 @@ function transform_to_model_space(p::Parameter{T,Exponential}, x::T) where T
     a + exp(c*(x-b))
 end
 
-transform_to_model_space(pvec::ParameterVector{T}, values::Vector{T}) where T = map(transform_to_model_space, pvec, values)
-transform_to_model_space(pvec::ParameterVector, values::Vector{S}) where S = map(transform_to_model_space, pvec, values)
+@inline function transform_to_model_space(pvec::ParameterVector{T}, values::Vector{T};
+                                          regime_switching::Bool = false) where T
+    if regime_switching
+        # Transform values in the first regime
+        output = similar(values)
+        plen   = length(pvec)
+        map!(transform_to_model_space, output, pvec, values[1:plen])
+
+        # Now transform values in the second regime and on
+        i = 0
+        for p in pvec
+            if !isempty(p.regimes)
+                for (k, v) in p.regimes[:value]
+                    if k != 1  # Skip the first regime.
+                        i += 1 # `values` stores regime values (after the first regime) beside each other.
+                        output[plen + i] = transform_to_model_space(p, values[plen + i])
+                    end
+                end
+            end
+        end
+
+        return output
+    else
+        return map(transform_to_model_space, pvec, values)
+    end
+end
+
+@inline function transform_to_model_space(pvec::ParameterVector, values::Vector{S};
+                                          regime_switching::Bool = false) where S
+
+    if regime_switching
+        # Transform values in the first regime
+        output = similar(values)
+        plen   = length(pvec)
+        map!(transform_to_model_space, output, pvec, values[1:plen])
+
+        # Now transform values in the second regime and on
+        i = 0
+        for p in pvec
+            if !isempty(p.regimes)
+                for (k, v) in p.regimes[:value]
+                    if k != 1  # Skip the first regime.
+                        i += 1 # `values` stores regime values (after the first regime) beside each other.
+                        output[plen + i] = transform_to_model_space(p, values[plen + i])
+                    end
+                end
+            end
+        end
+
+        return output
+    else
+        return map(transform_to_model_space, pvec, values)
+    end
+end
 
 """
 ```
@@ -837,8 +915,57 @@ function transform_to_real_line(p::Parameter{T,Exponential}, x::T = p.value) whe
     b + (1 ./ c) * log(x-a)
 end
 
-transform_to_real_line(pvec::ParameterVector{T}, values::Vector{T}) where T  = map(transform_to_real_line, pvec, values)
-transform_to_real_line(pvec::ParameterVector{T}) where T = map(transform_to_real_line, pvec)
+@inline function transform_to_real_line(pvec::ParameterVector{T}, values::Vector{T}; regime_switching::Bool = false) where T
+    if regime_switching
+        # Transform values in the first regime
+        output = similar(values)
+        plen   = length(pvec)
+        map!(transform_to_real_line, output, pvec, values[1:plen])
+
+        # Now transform values in the second regime and on
+        i = 0
+        for p in pvec
+            if !isempty(p.regimes)
+                for (k, v) in p.regimes[:value]
+                    if k != 1  # Skip the first regime.
+                        i += 1 # `values` stores regime values (after the first regime) beside each other.
+                        output[plen + i] = transform_to_real_line(p, values[plen + i])
+                    end
+                end
+            end
+        end
+
+        return output
+    else
+        map(transform_to_real_line, pvec, values)
+    end
+end
+@inline function transform_to_real_line(pvec::ParameterVector{T}; regime_switching::Bool = false) where T
+    if regime_switching
+        values = get_values(pvec) # regime-switching parameters returned by default
+
+        # Transform values in the first regime
+        plen = length(pvec) # since values is not passed, we can mutate values directly to avoid extra allocations
+        map!(transform_to_real_line, (@view values[1:plen]), pvec, values[1:plen])
+
+        # Now transform values in the second regime and on
+        i = 0
+        for p in pvec
+            if !isempty(p.regimes)
+                for (k, v) in p.regimes[:value]
+                    if k != 1  # Skip the first regime.
+                        i += 1 # `values` stores regime values (after the first regime) beside each other.
+                        values[plen + i] = transform_to_real_line(p, values[plen + i])
+                    end
+                end
+            end
+        end
+
+        return values
+    else
+        map(transform_to_real_line, pvec)
+    end
+end
 
 """
 ```

test/parameters.jl

@@ -39,7 +39,19 @@ tomodel_answers[3] = 1.
             @test differentiate_transform_to_real_line(u,u.value) != differentiate_transform_to_model_space(u,u.value)
         end
     end
+end
+
+@testset "Check type conversion for `value`, `valuebounds`, and `transform_parameterization`" begin
+    @info "The following warning is expected"
+    u1 = parameter(:σ_pist, 2.5230, (1, 5), (Float32(1), Float32(5)), Untransformed())
+    u2 = parameter(:σ_pist, 2.5230, (1., 5.), (1., 5.), Untransformed())
 
+    @test u1.value == u2.value
+    @test u1.valuebounds == u2.valuebounds
+    @test u1.transform_parameterization == u2.transform_parameterization
+    @test typeof(u1.value) == typeof(u2.value)
+    @test eltype(u1.valuebounds) == eltype(u2.valuebounds)
+    @test eltype(u1.transform_parameterization) == eltype(u2.transform_parameterization)
 end
 
 # probability

test/regimes.jl

@@ -1,6 +1,7 @@
 using Test, ModelConstructors
 
 # CURRENTLY ONLY TESTS VALUE, PRIOR, FIXED, AND VALUEBOUNDS SWITCHING, NO REGIME SWITCHING IN OTHER CASES
+# ALSO TESTS IF TRANSFORMS WORK CORRECTLY WITH REGIME-SWITCHING
 
 @info "The following error 'get_regime_val(), Input Error: No regime 3' is expected."
 @testset "Regime switching with parameters" begin
@@ -67,6 +68,24 @@ using Test, ModelConstructors
     ModelConstructors.set_regime_fixed!(u, 2, false; update_valuebounds = (10., 11.))
     @test u.regimes[:valuebounds][1] == (10., 11.)
     @test u.regimes[:valuebounds][2] == (10., 11.)
+
+    # test transform_to_real_line
+    ModelConstructors.toggle_regime!(uvec, 1)
+    set_regime_val!(uvec[1], 2, 1.; override_bounds = true) # break the valuebounds for now
+    set_regime_val!(uvec[2], 2, 1.; override_bounds = true)
+    values = ModelConstructors.get_values(uvec)
+    real_vals_true = similar(values)
+    real_vals_true[1] = ModelConstructors.transform_to_real_line(uvec[1], uvec[1].regimes[:value][1])
+    real_vals_true[2] = ModelConstructors.transform_to_real_line(uvec[2], uvec[2].regimes[:value][1])
+    real_vals_true[3] = ModelConstructors.transform_to_real_line(uvec[1], uvec[1].regimes[:value][2])
+    real_vals_true[4] = ModelConstructors.transform_to_real_line(uvec[2], uvec[2].regimes[:value][2])
+    real_vals1 = transform_to_real_line(uvec; regime_switching = true)
+    real_vals2 = transform_to_real_line(uvec, values; regime_switching = true)
+    @test real_vals1 == real_vals2 == real_vals_true
+
+    # test transform_to_model_space
+    model_vals = transform_to_model_space(uvec, real_vals1; regime_switching = true)
+    @test model_vals == values
 end
 
 @testset "Regime switching with parameters when model regimes are different" begin