Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

CompatHelper: bump compat for "StatsModels" to "0.7" #62

Merged
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
14 changes: 7 additions & 7 deletions Project.toml
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
name = "JellyMe4"
uuid = "19ac8677-9a15-4623-9afd-84acc6165ce7"
authors = ["Phillip Alday <[email protected]>"]
version = "0.2.7"
version = "0.2.8"

[deps]
CategoricalArrays = "324d7699-5711-5eae-9e2f-1d82baa6b597"
Expand All @@ -17,16 +17,16 @@ StatsModels = "3eaba693-59b7-5ba5-a881-562e759f1c8d"
Tables = "bd369af6-aec1-5ad0-b16a-f7cc5008161c"

[compat]
CategoricalArrays = "0.7, 0.8, 0.9, 0.10"
DataFrames = "0.20, 0.21, 0.22, 1.0"
Distributions = "0.20, 0.22, 0.23, 0.24, 0.25"
CategoricalArrays = "0.10"
DataFrames = "1"
Distributions = "0.25"
GLM = "1.3"
MixedModels = "3.0, 4.0"
MixedModels = "4"
RCall = "0.13"
RegressionFormulae = "0.1"
StatsBase = "0.31, 0.32, 0.33"
StatsModels = "0.6.8"
Tables = "0.2, 1.0"
StatsModels = "0.6, 0.7"
Tables = "1"
julia = "1.6"

[extras]
Expand Down
2 changes: 1 addition & 1 deletion deps/build.jl
Original file line number Diff line number Diff line change
@@ -1 +1 @@
# should we try to install lme4 or even afex automatically at build time?
# should we try to install lme4 or even afex automatically at build time?
2 changes: 1 addition & 1 deletion src/glmerMod.jl
Original file line number Diff line number Diff line change
Expand Up @@ -119,7 +119,7 @@ function rcopy(::Type{GeneralizedLinearMixedModel}, s::Ptr{S4Sxp})

if length(θ) != length(m.θ)
@error """You're probably using || in R with a categorical variable,
whose translation is currently unsupported with MixedModels 3.0."""
whose translation is currently unsupported with recent MixedModels releases"""
throw(ArgumentError("Parameter vectors in R and Julia are different sizes."))
end

Expand Down
2 changes: 1 addition & 1 deletion src/lmerMod.jl
Original file line number Diff line number Diff line change
Expand Up @@ -57,7 +57,7 @@ function rcopy(::Type{LinearMixedModel}, s::Ptr{S4Sxp})

if length(θ) != length(m.θ)
@error """You're probably using || in R with a categorical variable,
whose translation is currently unsupported with MixedModels 3.0."""
whose translation is currently unsupported with recent MixedModels releases"""
throw(ArgumentError("Parameter vectors in R and Julia are different sizes."))
end

Expand Down
22 changes: 18 additions & 4 deletions src/utilities.jl
Original file line number Diff line number Diff line change
Expand Up @@ -41,16 +41,30 @@ function get_r_contrasts(rdf)
data = rcopy(rdf)
# get categorical columns
cnames = [c for c in propertynames(data) if typeof(data[!, c]) <: CategoricalArray]
return Dict(c => HypothesisCoding(pinv(rcopy(R"contrasts($(rdf[c]))"));
labels=rcopyarray(R"colnames(contrasts($(rdf[c])))"))
for c in cnames)
pairs = []
for c in cnames
@debug "" c
levels = rcopyarray(R"rownames(contrasts($(rdf[c])))")
@debug "" levels
if levels == [nothing]
# for some contrasts, the rownames aren't defined
# we don't prefer this method because it's possible
# that it's not quite right when people are using really
# fancy contrasts and it can fail
levels = rcopyarray(R"levels($(rdf[c]))")
end
labels = rcopyarray(R"colnames(contrasts($(rdf[c])))")
@debug "" labels
push!(pairs,
c => HypothesisCoding(pinv(rcopy(R"contrasts($(rdf[c]))")); levels, labels))
end
return Dict(pairs)
end

"""
set the contrasts on an R dataframe
this of course has a side effect: it changes the R dataframe
"""

function set_r_contrasts!(rdfname, formula)
fixefform = first(formula.rhs)

Expand Down
12 changes: 6 additions & 6 deletions test/glmerMod.jl
Original file line number Diff line number Diff line change
Expand Up @@ -132,7 +132,7 @@ logistic(x) = 1 / (1 + exp(-x))
# we max out at 1 scalar RE for the nAGQ tests
jlmm = GLMM(@formula(r2 ~ 1 + anger + gender + btype + situ + (1 | subj)),
dat, Bernoulli())
fit!(jlmm; fast=false, nAGQ=9)
fit!(jlmm; fast=false, nAGQ=9, progress=false)
jm = (jlmm, dat)
@rput jm
# @test_warn Regex(".*categorical.*") @rput jm;
Expand All @@ -147,7 +147,7 @@ logistic(x) = 1 / (1 + exp(-x))
# we max out at 1 scalar RE for the nAGQ tests
jlmm = GLMM(@formula(r2 ~ 1 + anger + gender + btype + situ + (1 | subj)),
dat, Bernoulli())
fit!(jlmm)
fit!(jlmm; progress=false)
jm = Tuple([jlmm, dat])
@rput jm

Expand All @@ -159,7 +159,7 @@ logistic(x) = 1 / (1 + exp(-x))
@testset "columntable" begin
jlmm = GLMM(@formula(r2 ~ 1 + anger + gender + btype + situ + (1 | subj)),
dat, Bernoulli())
fit!(jlmm)
fit!(jlmm; progress=false)
jm = Tuple([jlmm, columntable(dat)])
@rput jm
end
Expand All @@ -170,7 +170,7 @@ logistic(x) = 1 / (1 + exp(-x))
dat = dataset(:cbpp)
dat.rate = dat.incid ./ dat.hsz
jlmm = fit(MixedModel, @formula(rate ~ 1 + period + (1 | herd)),
dat, Binomial(); wts=float(dat.hsz), fast=true)
dat, Binomial(); wts=float(dat.hsz), fast=true, progress=false)

jm = (jlmm, dat)
@rput jm
Expand All @@ -192,7 +192,7 @@ logistic(x) = 1 / (1 + exp(-x))
dat, Poisson())

jm = (jlmm, dat)
fit!(jlmm; fast=true) # problems with this one in fast=false
fit!(jlmm; fast=true, progress=false) # problems with this one in fast=false
@rput jm
# @test_warn Regex(".*categorical.*") @rput jm;
@test rcopy(R"""jm@devcomp$dims["nAGQ"]""") == 0
Expand Down Expand Up @@ -228,7 +228,7 @@ logistic(x) = 1 / (1 + exp(-x))
jlmm = fit(MixedModel,
@formula(r2 ~ 1 + anger + gender + btype + situ + (1 | subj) +
(1 | item)),
dat, Bernoulli();
dat, Bernoulli(); progress=false,
contrasts=Dict(:gender => EffectsCoding(),
:btypes => EffectsCoding()))
jm = (jlmm, dat)
Expand Down
44 changes: 23 additions & 21 deletions test/lmerMod.jl
Original file line number Diff line number Diff line change
@@ -1,9 +1,10 @@
using RCall, MixedModels, Test
using RCall, JellyMe4, MixedModels, Test
using StatsBase: zscore
using StatsModels: SeqDiffCoding
using Tables: columntable

import JellyMe4: _set_lmer, _set_afex_installed
using JellyMe4: _set_lmer, _set_afex_installed
using MixedModels: dataset

const LMM = LinearMixedModel
const GLMM = GeneralizedLinearMixedModel
Expand All @@ -25,22 +26,22 @@ const GLMM = GeneralizedLinearMixedModel

# this is available in MixedModels.dataset(:sleepstudy) but with different
# capitalization than in R
sleepstudy = rcopy(R"sleepstudy")
sleepstudy = rcopy(R"lme4::sleepstudy")

kb07 = dataset(:kb07)

@testset "get lmerMod" begin
### from R ###

jlmm = fit!(LMM(@formula(Reaction ~ 1 + Days + (1 + Days | Subject)), sleepstudy);
REML=false)
rlmm = rcopy(R"m <- lme4::lmer(Reaction ~ 1 + Days + (1 + Days|Subject),sleepstudy,REML=FALSE)")
REML=false, progress=false)
rlmm = rcopy(R"m <- lme4::lmer(Reaction ~ 1 + Days + (1 + Days|Subject),lme4::sleepstudy,REML=FALSE)")

@test jlmm.θ ≈ rlmm.θ atol = 0.001
@test objective(jlmm) ≈ objective(rlmm) atol = 0.001
@test fixef(jlmm) ≈ fixef(rlmm) atol = 0.001

jlmm = refit!(jlmm; REML=true)
jlmm = refit!(jlmm; REML=true, progress=false)
rlmm = rcopy(R"update(m, REML=TRUE)")

@test jlmm.θ ≈ rlmm.θ atol = 0.001
Expand All @@ -49,14 +50,14 @@ const GLMM = GeneralizedLinearMixedModel

@testset "merModLmerTest" begin
jlmm = fit!(LMM(@formula(Reaction ~ 1 + Days + (1 + Days | Subject)),
sleepstudy); REML=false)
sleepstudy); REML=false, progress=false)
rlmm = rcopy(R"m <- lmerTest::lmer(Reaction ~ 1 + Days + (1 + Days|Subject),sleepstudy,REML=FALSE)")

@test jlmm.θ ≈ rlmm.θ atol = 0.001
@test objective(jlmm) ≈ objective(rlmm) atol = 0.001
@test fixef(jlmm) ≈ fixef(rlmm) atol = 0.001

jlmm = refit!(jlmm; REML=true)
jlmm = refit!(jlmm; REML=true, progress=false)
rlmm = rcopy(R"update(m, REML=TRUE)")

@test jlmm.θ ≈ rlmm.θ atol = 0.001
Expand All @@ -70,7 +71,7 @@ const GLMM = GeneralizedLinearMixedModel
# the scalar-vector distinction for θ is missing in R
@testset "scalar RE" begin
jlmm = fit!(LMM(@formula(Reaction ~ 1 + Days + (1 | Subject)), sleepstudy);
REML=false)
REML=false, progress=false)
rlmm = rcopy(R"m <- lme4::lmer(Reaction ~ 1 + Days + (1|Subject),sleepstudy,REML=FALSE)")

@test jlmm.θ ≈ rlmm.θ atol = 0.001
Expand All @@ -91,7 +92,7 @@ const GLMM = GeneralizedLinearMixedModel

@testset "contrasts" begin
reval("""
data(cake)
cake <- lme4::cake
cake\$rr <- with(cake, replicate:recipe)
""")
rlmm = rcopy(R"fm1 <- lme4::lmer(angle ~ recipe * temperature + (1|rr), cake, REML= FALSE)")
Expand All @@ -103,7 +104,7 @@ const GLMM = GeneralizedLinearMixedModel

@testset "caret" begin
reval("""
data(cake)
cake <- lme4::cake
cake\$rr <- with(cake, replicate:recipe)
""")
rlmm = rcopy(R"fm1 <- lme4::lmer(angle ~ (recipe + temperature)^2 + (1|rr), cake, REML= FALSE)")
Expand All @@ -130,7 +131,7 @@ const GLMM = GeneralizedLinearMixedModel

rlmm = rcopy(R"m <- lme4::lmer(Reaction ~ 1 + Days + (1 + Days||Subject),sleepstudy,REML=FALSE)")
jlmm = fit!(LMM(@formula(Reaction ~ 1 + Days + zerocorr(1 + Days | Subject)),
sleepstudy); REML=false)
sleepstudy); REML=false, progress=false)
# as a cheat for comparing the covariance matrices, we use PCA
@test only(rlmm.rePCA) ≈ only(jlmm.rePCA) atol = 0.05
end
Expand All @@ -151,7 +152,7 @@ const GLMM = GeneralizedLinearMixedModel
@test rcopy(R"fitted(jm)") ≈ fitted(jlmm)
@test rcopy(R"REMLcrit(jm)") ≈ objective(jlmm)

refit!(jlmm; REML=false)
refit!(jlmm; REML=false, progress=false)
jm = (jlmm, sleepstudy)
@rput jm
@test rcopy(R"fitted(jm)") ≈ fitted(jlmm)
Expand All @@ -168,14 +169,14 @@ const GLMM = GeneralizedLinearMixedModel
R"m <- lme4::lmer(log10(Reaction) ~ 1 + log(Days2) + (1 + log(Days2)|Subject),sleepstudy,REML=FALSE)"
jlmm = fit!(LMM(@formula(log10(Reaction) ~ 1 + log(Days2) +
(1 + log(Days2) | Subject)),
sleepstudy); REML=false)
sleepstudy); REML=false, progress=false)
jm = (jlmm, sleepstudy)
@rput jm
@test rcopy(R"fitted(jm)") ≈ fitted(jlmm)
@test rcopy(R"deviance(jm)") ≈ objective(jlmm)

jlmm = fit!(LMM(@formula(Reaction ~ 1 + round(Days) + (1 | Subject)),
sleepstudy); REML=false)
sleepstudy); REML=false, progress=false)
jm = (jlmm, sleepstudy)

@test_throws ArgumentError (@rput jm)
Expand All @@ -184,7 +185,7 @@ const GLMM = GeneralizedLinearMixedModel
@testset "sorting by n BLUPs vs. n groups" begin
jlmm = fit(LMM,
@formula(rt_trunc ~ 1 + (1 | subj) + (1 + load + prec + spkr | item)),
kb07)
kb07; progress=false)
jm = (jlmm, kb07)
@rput jm
@test rcopy(R"fitted(jm)") ≈ fitted(jlmm)
Expand All @@ -202,7 +203,8 @@ const GLMM = GeneralizedLinearMixedModel
cake
"""))
jlmm = fit(MixedModel, @formula(angle ~ recipe * temperature + (1 | rr)),
cake; REML=false, contrasts=Dict(:temperature => SeqDiffCoding()))
cake; REML=false, progress=false,
contrasts=Dict(:temperature => SeqDiffCoding()))
jm = (jlmm, cake)
@rput jm
@test fixef(jlmm) ≈ rcopy(R"fixef(jm)")
Expand All @@ -212,7 +214,7 @@ const GLMM = GeneralizedLinearMixedModel
machines = rcopy(R"as.data.frame(nlme::Machines)")
jlmm = fit(MixedModel,
@formula(score ~ 1 + Machine + (1 + fulldummy(Machine) | Worker)),
machines)
machines; progress=false)
rlmm = (jlmm, machines)
@rput rlmm
rlmmrepca = rcopy(R"summary(rePCA(rlmm))$Worker$importance[3,]")
Expand All @@ -227,7 +229,7 @@ const GLMM = GeneralizedLinearMixedModel
_set_afex_installed(false)

jlmm = fit!(LMM(@formula(Reaction ~ 1 + Days + zerocorr(1 + Days | Subject)),
sleepstudy); REML=false)
sleepstudy); REML=false, progress=false)
rlmm = (jlmm, sleepstudy)
@rput rlmm
@test rcopy(R"""!is(rlmm,"merModLmerTest")""")
Expand All @@ -241,7 +243,7 @@ const GLMM = GeneralizedLinearMixedModel

jlmm = fit(MixedModel,
@formula(score ~ 1 + Machine + zerocorr(0 + Machine | Worker)),
machines)
machines; progress=false)

@testset "afex pre-enabled" begin
_set_lmer("afex::lmer_alt")
Expand Down Expand Up @@ -279,7 +281,7 @@ const GLMM = GeneralizedLinearMixedModel
slp[!, :obs] .= 1:nrow(slp)
fm1 = fit(MixedModel,
@formula(reaction ~ 1 + days + (1 | obs)),
slp;
slp; progress=false,
contrasts=Dict(:obs => Grouping()))
rfm1 = (fm1, slp)
warning = r"""Warning: number of levels of each grouping factor must be < number of observations \(problems: obs\)
Expand Down
4 changes: 2 additions & 2 deletions test/merMod.jl
Original file line number Diff line number Diff line change
Expand Up @@ -6,9 +6,9 @@ const GLMM = GeneralizedLinearMixedModel
@testset "merMod" begin
# this is available in MixedModels.dataset(:sleepstudy) but with different
# capitalization than in R
sleepstudy = rcopy(R"sleepstudy")
sleepstudy = rcopy(R"lme4::sleepstudy")
jlmm = fit!(LMM(@formula(Reaction ~ 1 + round(Days) + (1 | Subject)), sleepstudy);
REML=false)
REML=false, progress=false)
@testset "bare model" begin
@test_throws ArgumentError (@rput jlmm)
end
Expand Down