rename n_obs() -> nsamples()

alyst · alyst · commit 273cb02455e5 · 2024-06-16T17:09:23.000-07:00
diff --git a/src/StructuralEquationModels.jl b/src/StructuralEquationModels.jl
@@ -80,7 +80,6 @@ include("frontend/fit/fitmeasures/BIC.jl")
 include("frontend/fit/fitmeasures/chi2.jl")
 include("frontend/fit/fitmeasures/df.jl")
 include("frontend/fit/fitmeasures/minus2ll.jl")
-include("frontend/fit/fitmeasures/n_obs.jl")
 include("frontend/fit/fitmeasures/p.jl")
 include("frontend/fit/fitmeasures/RMSEA.jl")
 include("frontend/fit/fitmeasures/n_man.jl")
@@ -167,7 +166,7 @@ export AbstractSem,
     df,
     fit_measures,
     minus2ll,
-    n_obs,
+    nsamples,
     p_value,
     RMSEA,
     n_man,
diff --git a/src/frontend/common.jl b/src/frontend/common.jl
@@ -45,4 +45,13 @@ parind[:param_name]
 
 See also [`params`](@ref).
 """
-param_indices(semobj) = Dict(par => i for (i, par) in enumerate(params(semobj)))
+param_indices(semobj) = Dict(par => i for (i, par) in enumerate(params(semobj)))
+
+"""
+    nsamples(semobj)
+
+Return the number of samples (observed data points).
+
+For ensemble models, return the sum over all submodels.
+"""
+function nsamples end
diff --git a/src/frontend/fit/SemFit.jl b/src/frontend/fit/SemFit.jl
@@ -50,6 +50,7 @@ end
 
 params(fit::SemFit) = params(fit.model)
 nparams(fit::SemFit) = nparams(fit.model)
+nsamples(fit::SemFit) = nsamples(fit.model)
 
 # access fields
 minimum(sem_fit::SemFit) = sem_fit.minimum
diff --git a/src/frontend/fit/fitmeasures/BIC.jl b/src/frontend/fit/fitmeasures/BIC.jl
@@ -3,4 +3,4 @@
 
 Return the bayesian information criterion.
 """
-BIC(sem_fit::SemFit) = minus2ll(sem_fit) + log(n_obs(sem_fit)) * nparams(sem_fit)
+BIC(sem_fit::SemFit) = minus2ll(sem_fit) + log(nsamples(sem_fit)) * nparams(sem_fit)
diff --git a/src/frontend/fit/fitmeasures/RMSEA.jl b/src/frontend/fit/fitmeasures/RMSEA.jl
@@ -8,13 +8,13 @@ function RMSEA end
 RMSEA(sem_fit::SemFit) = RMSEA(sem_fit, sem_fit.model)
 
 RMSEA(sem_fit::SemFit, model::AbstractSemSingle) =
-    RMSEA(df(sem_fit), χ²(sem_fit), n_obs(sem_fit))
+    RMSEA(df(sem_fit), χ²(sem_fit), nsamples(sem_fit))
 
 RMSEA(sem_fit::SemFit, model::SemEnsemble) =
-    sqrt(length(model.sems)) * RMSEA(df(sem_fit), χ²(sem_fit), n_obs(sem_fit))
+    sqrt(length(model.sems)) * RMSEA(df(sem_fit), χ²(sem_fit), nsamples(sem_fit))
 
-function RMSEA(df, chi2, n_obs)
-    rmsea = (chi2 - df) / (n_obs * df)
+function RMSEA(df, chi2, nsamples)
+    rmsea = (chi2 - df) / (nsamples * df)
     rmsea > 0 ? nothing : rmsea = 0
     return sqrt(rmsea)
 end
diff --git a/src/frontend/fit/fitmeasures/chi2.jl b/src/frontend/fit/fitmeasures/chi2.jl
@@ -14,11 +14,11 @@ Return the χ² value.
 
 # RAM + SemML
 χ²(lossfun::SemML, fit::SemFit, model::AbstractSemSingle) =
-    (n_obs(fit) - 1) *
+    (nsamples(fit) - 1) *
     (fit.minimum - logdet(obs_cov(observed(model))) - n_man(observed(model)))
 
 # bollen, p. 115, only correct for GLS weight matrix
-χ²(lossfun::SemWLS, fit::SemFit, model::AbstractSemSingle) = (n_obs(fit) - 1) * fit.minimum
+χ²(lossfun::SemWLS, fit::SemFit, model::AbstractSemSingle) = (nsamples(fit) - 1) * fit.minimum
 
 # FIML
 function χ²(lossfun::SemFIML, fit::SemFit, model::AbstractSemSingle)
@@ -51,15 +51,15 @@ function χ²(fit::SemFit, models::SemEnsemble)
 end
 
 function χ²(lossfun::SemWLS, fit::SemFit, models::SemEnsemble)
-    return (sum(n_obs, models.sems) - 1) * fit.minimum
+    return (nsamples(models) - 1) * fit.minimum
 end
 
 function χ²(lossfun::SemML, fit::SemFit, models::SemEnsemble)
     G = sum(zip(models.weights, models.sems)) do (w, model)
         data = observed(model)
         w * (logdet(obs_cov(data)) + n_man(data))
     end
-    return (sum(n_obs, models.sems) - 1) * (fit.minimum - G)
+    return (nsamples(models) - 1) * (fit.minimum - G)
 end
 
 function χ²(lossfun::SemFIML, fit::SemFit, models::SemEnsemble)
diff --git a/src/frontend/fit/fitmeasures/minus2ll.jl b/src/frontend/fit/fitmeasures/minus2ll.jl
@@ -22,7 +22,7 @@ minus2ll(minimum::Number, model::AbstractSemSingle) =
 # SemML ------------------------------------------------------------------------------------
 function minus2ll(lossfun::SemML, minimum::Number, model::AbstractSemSingle)
     obs = observed(model)
-    return n_obs(obs) * (minimum + log(2π) * n_man(obs))
+    return nsamples(obs) * (minimum + log(2π) * n_man(obs))
 end
 
 # WLS --------------------------------------------------------------------------------------
@@ -32,8 +32,8 @@ minus2ll(lossfun::SemWLS, minimum::Number, model::AbstractSemSingle) = missing
 # -2ll = (∑ log(2π)*(nᵢ + mᵢ)) + F*n
 function minus2ll(lossfun::SemFIML, minimum::Number, model::AbstractSemSingle)
     obs = observed(model)::SemObservedMissing
-    F = minimum * n_obs(obs)
-    F += log(2π) * sum(pat -> n_obs(pat) * nobserved_vars(pat), obs.patterns)
+    F = minimum * nsamples(obs)
+    F += log(2π) * sum(pat -> nsamples(pat) * nobserved_vars(pat), obs.patterns)
     return F
 end
 
@@ -44,7 +44,7 @@ function minus2ll(observed::SemObservedMissing)
 
     F = 0.0
     for pat in observed.patterns
-        nᵢ = n_obs(pat)
+        nᵢ = nsamples(pat)
         # implied covariance/mean
         Σᵢ = Symmetric(Σ[pat.obs_mask, pat.obs_mask])
 
@@ -53,13 +53,13 @@ function minus2ll(observed::SemObservedMissing)
         μ_diffᵢ = pat.obs_mean - μ[pat.obs_mask]
 
         F_pat = ld + dot(μ_diffᵢ, Σᵢ⁻¹, μ_diffᵢ)
-        if n_obs(pat) > 1
+        if nsamples(pat) > 1
             F_pat += dot(pat.obs_cov, Σᵢ⁻¹)
         end
-        F += (F_pat + log(2π) * nobserved_vars(pat)) * n_obs(pat)
+        F += (F_pat + log(2π) * nobserved_vars(pat)) * nsamples(pat)
     end
 
-    #F *= n_obs(observed)
+    #F *= nsamples(observed)
     return F
 end
 
diff --git a/src/frontend/fit/fitmeasures/n_obs.jl b/src/frontend/fit/fitmeasures/n_obs.jl
diff --git a/src/frontend/fit/standard_errors/hessian.jl b/src/frontend/fit/standard_errors/hessian.jl
@@ -42,13 +42,13 @@ function H_scaling(model::AbstractSemSingle)
     return H_scaling(model.loss.functions[1], model)
 end
 
-H_scaling(lossfun::SemML, model::AbstractSemSingle) = 2 / (n_obs(model) - 1)
+H_scaling(lossfun::SemML, model::AbstractSemSingle) = 2 / (nsamples(model) - 1)
 
 function H_scaling(lossfun::SemWLS, model::AbstractSemSingle)
     @warn "Standard errors for WLS are only correct if a GLS weight matrix (the default) is used."
-    return 2 / (n_obs(model) - 1)
+    return 2 / (nsamples(model) - 1)
 end
 
-H_scaling(lossfun::SemFIML, model::AbstractSemSingle) = 2 / n_obs(model)
+H_scaling(lossfun::SemFIML, model::AbstractSemSingle) = 2 / nsamples(model)
 
-H_scaling(model::SemEnsemble) = 2 / n_obs(model)
+H_scaling(model::SemEnsemble) = 2 / nsamples(model)
diff --git a/src/frontend/fit/summary.jl b/src/frontend/fit/summary.jl
@@ -17,7 +17,7 @@ function sem_summary(
     println("No. iterations/evaluations:  $(n_iterations(sem_fit))")
     print("\n")
     println("Number of parameters:        $(nparams(sem_fit))")
-    println("Number of observations:      $(n_obs(sem_fit))")
+    println("Number of data samples:      $(nsamples(sem_fit))")
     print("\n")
     printstyled(
         "----------------------------------- Model ----------------------------------- \n";
diff --git a/src/frontend/specification/Sem.jl b/src/frontend/specification/Sem.jl
@@ -28,6 +28,11 @@ vars(sem::AbstractSemSingle) = vars(sem.imply.ram_matrices)
 observed_vars(sem::AbstractSemSingle) = observed_vars(sem.imply.ram_matrices)
 latent_vars(sem::AbstractSemSingle) = latent_vars(sem.imply.ram_matrices)
 
+nsamples(sem::AbstractSemSingle) = nsamples(sem.observed)
+
+# sum of samples in all sub-models
+nsamples(ensemble::SemEnsemble) = sum(nsamples, ensemble.sems)
+
 function SemFiniteDiff(;
     observed::O = SemObservedData,
     imply::I = RAM,
diff --git a/src/loss/ML/FIML.jl b/src/loss/ML/FIML.jl
@@ -41,25 +41,25 @@ end
 
 function objective(fiml::SemFIMLPattern{T}, pat::SemObservedMissingPattern) where {T}
     F = fiml.logdet[] + dot(fiml.μ_diff, fiml.Σ⁻¹, fiml.μ_diff)
-    if n_obs(pat) > 1
+    if nsamples(pat) > 1
         F += dot(pat.obs_cov, fiml.Σ⁻¹)
-        F *= n_obs(pat)
+        F *= nsamples(pat)
     end
     return F
 end
 
 function gradient!(JΣ, Jμ, fiml::SemFIMLPattern, pat::SemObservedMissingPattern)
     Σ⁻¹ = Symmetric(fiml.Σ⁻¹)
     μ_diff⨉Σ⁻¹ = fiml.μ_diff' * Σ⁻¹
-    if n_obs(pat) > 1
+    if nsamples(pat) > 1
         JΣ_pat = Σ⁻¹ * (I - pat.obs_cov * Σ⁻¹ - fiml.μ_diff * μ_diff⨉Σ⁻¹)
-        JΣ_pat .*= n_obs(pat)
+        JΣ_pat .*= nsamples(pat)
     else
         JΣ_pat = Σ⁻¹ * (I - fiml.μ_diff * μ_diff⨉Σ⁻¹)
     end
     @inbounds vec(JΣ)[fiml.∇ind] .+= vec(JΣ_pat)
 
-    lmul!(2 * n_obs(pat), μ_diff⨉Σ⁻¹)
+    lmul!(2 * nsamples(pat), μ_diff⨉Σ⁻¹)
     @inbounds Jμ[pat.obs_mask] .+= μ_diff⨉Σ⁻¹'
     return nothing
 end
@@ -133,7 +133,7 @@ function evaluate!(
 
     prepare!(fiml, model)
 
-    scale = inv(n_obs(observed(model)))
+    scale = inv(nsamples(observed(model)))
     isnothing(objective) ||
         (objective = scale * F_FIML(eltype(params), fiml, observed(model), model))
     isnothing(gradient) ||
diff --git a/src/observed/EM.jl b/src/observed/EM.jl
@@ -41,7 +41,7 @@ function em_mvn(
         if nmissed_vars(pat) == 0
             𝔼x_full .+= sum(pat.data, dims = 2)
             mul!(𝔼xxᵀ_full, pat.data, pat.data', 1, 1)
-            nobs_full += n_obs(pat)
+            nobs_full += nsamples(pat)
         end
     end
     if nobs_full == 0
@@ -124,10 +124,10 @@ function em_step!(
         μ₀o = μ₀[o]
 
         # get pattern observations
-        nobs = !isnothing(max_nobs_em) ? min(max_nobs_em, n_obs(pat)) : n_obs(pat)
+        nobs = !isnothing(max_nobs_em) ? min(max_nobs_em, nsamples(pat)) : nsamples(pat)
         zo =
-            nobs < n_obs(pat) ?
-            pat.data[:, sort!(sample(1:n_obs(pat), nobs, replace = false))] : copy(pat.data)
+            nobs < nsamples(pat) ?
+            pat.data[:, sort!(sample(1:nsamples(pat), nobs, replace = false))] : copy(pat.data)
         zo .-= μ₀o # subtract current mean from observations
 
         𝔼zo = sum(zo, dims = 2)
@@ -202,7 +202,7 @@ end
 # use μ and Σ of full cases
 function start_em_observed(patterns::AbstractVector{<:SemObservedMissingPattern}; kwargs...)
     fullpat = patterns[1]
-    if (nmissed_vars(fullpat) == 0) && (n_obs(fullpat) > 1)
+    if (nmissed_vars(fullpat) == 0) && (nsamples(fullpat) > 1)
         μ = copy(fullpat.obs_mean)
         Σ = copy(parent(fullpat.obs_cov))
         if !isposdef(Σ)
diff --git a/src/observed/covariance.jl b/src/observed/covariance.jl
@@ -10,7 +10,7 @@ function SemObservedCovariance(;
     spec_colnames::Union{AbstractVector{Symbol}, Nothing} = nothing,
     obs_mean::Union{AbstractVector, Nothing} = nothing,
     meanstructure::Bool = false,
-    n_obs::Integer,
+    nsamples::Integer,
     kwargs...,
 )
     if !meanstructure && !isnothing(obs_mean)
@@ -33,5 +33,5 @@ function SemObservedCovariance(;
         isnothing(obs_mean) || (obs_mean = obs_mean[obs2spec_perm])
     end
 
-    return SemObservedData(nothing, Symmetric(obs_cov), obs_mean, size(obs_cov, 1), n_obs)
+    return SemObservedData(nothing, Symmetric(obs_cov), obs_mean, size(obs_cov, 1), nsamples)
 end
diff --git a/src/observed/data.jl b/src/observed/data.jl
@@ -18,7 +18,7 @@ For observed data without missings.
 
 # Extended help
 ## Interfaces
-- `n_obs(::SemObservedData)` -> number of observed data points
+- `nsamples(::SemObservedData)` -> number of observed data points
 - `n_man(::SemObservedData)` -> number of manifest variables
 
 - `get_data(::SemObservedData)` -> observed data
@@ -42,7 +42,7 @@ struct SemObservedData{A, B, C} <: SemObserved
     obs_cov::B
     obs_mean::C
     n_man::Int
-    n_obs::Int
+    nsamples::Int
 end
 
 # error checks
@@ -112,7 +112,7 @@ end
 ### Recommended methods
 ############################################################################################
 
-n_obs(observed::SemObservedData) = observed.n_obs
+nsamples(observed::SemObservedData) = observed.nsamples
 n_man(observed::SemObservedData) = observed.n_man
 
 ############################################################################################
diff --git a/src/observed/missing.jl b/src/observed/missing.jl
@@ -20,7 +20,7 @@ For observed data with missing values.
 
 # Extended help
 ## Interfaces
-- `n_obs(::SemObservedMissing)` -> number of observed data points
+- `nsamples(::SemObservedMissing)` -> number of observed data points
 - `n_man(::SemObservedMissing)` -> number of manifest variables
 
 - `get_data(::SemObservedMissing)` -> observed data
@@ -41,7 +41,7 @@ struct SemObservedMissing{A <: AbstractMatrix, P <: SemObservedMissingPattern, T
        SemObserved
     data::A
     n_man::Int
-    n_obs::Int
+    nsamples::Int
     patterns::Vector{P}
 
     obs_cov::Matrix{T}
@@ -98,7 +98,7 @@ function SemObservedMissing(;
         data = Matrix(data)
     end
 
-    n_obs, n_man = size(data)
+    nsamples, n_man = size(data)
 
     # detect all different missing patterns with their row indices
     pattern_to_rows = Dict{BitVector, Vector{Int}}()
@@ -117,10 +117,10 @@ function SemObservedMissing(;
 
     em_cov, em_mean = em_mvn(patterns; kwargs...)
 
-    return SemObservedMissing(data, n_man, n_obs, patterns, em_cov, em_mean)
+    return SemObservedMissing(data, n_man, nsamples, patterns, em_cov, em_mean)
 end
 
-n_obs(observed::SemObservedMissing) = observed.n_obs
+nsamples(observed::SemObservedMissing) = observed.nsamples
 n_man(observed::SemObservedMissing) = observed.n_man
 
 obs_cov(observed::SemObservedMissing) = observed.obs_cov
diff --git a/src/observed/missing_pattern.jl b/src/observed/missing_pattern.jl
@@ -42,7 +42,7 @@ function SemObservedMissingPattern(
 end
 
 n_man(pat::SemObservedMissingPattern) = length(pat.obs_mask)
-n_obs(pat::SemObservedMissingPattern) = length(pat.rows)
+nsamples(pat::SemObservedMissingPattern) = length(pat.rows)
 
 nobserved_vars(pat::SemObservedMissingPattern) = pat.nobserved
 nmissed_vars(pat::SemObservedMissingPattern) = pat.nmissed
diff --git a/src/types.jl b/src/types.jl
@@ -228,8 +228,8 @@ function SemEnsemble(models...; optimizer = SemOptimizerOptim, weights = nothing
     # default weights
 
     if isnothing(weights)
-        nobs_total = sum(n_obs, models)
-        weights = [n_obs(model) / nobs_total for model in models]
+        nsamples_total = sum(nsamples, models)
+        weights = [nsamples(model) / nsamples_total for model in models]
     end
 
     # check parameters equality
diff --git a/test/examples/political_democracy/by_parts.jl b/test/examples/political_democracy/by_parts.jl
@@ -41,7 +41,7 @@ model_ridge = Sem(observed, imply_ram, SemLoss(ml, ridge), optimizer_obj)
 model_constant = Sem(observed, imply_ram, SemLoss(ml, constant), optimizer_obj)
 
 model_ml_weighted =
-    Sem(observed, imply_ram, SemLoss(ml; loss_weights = [n_obs(model_ml)]), optimizer_obj)
+    Sem(observed, imply_ram, SemLoss(ml; loss_weights = [nsamples(model_ml)]), optimizer_obj)
 
 ############################################################################################
 ### test gradients
@@ -101,7 +101,7 @@ end
     solution_ml = sem_fit(model_ml)
     solution_ml_weighted = sem_fit(model_ml_weighted)
     @test solution(solution_ml) ≈ solution(solution_ml_weighted) rtol = 1e-3
-    @test n_obs(model_ml) * StructuralEquationModels.minimum(solution_ml) ≈
+    @test nsamples(model_ml) * StructuralEquationModels.minimum(solution_ml) ≈
           StructuralEquationModels.minimum(solution_ml_weighted) rtol = 1e-6
 end
 
diff --git a/test/examples/political_democracy/constructor.jl b/test/examples/political_democracy/constructor.jl
diff --git a/test/unit_tests/data_input_formats.jl b/test/unit_tests/data_input_formats.jl
