Fit a (mixed-effects) regression model in Julia
Usage
jlm(formula, data, family = "gaussian", contrasts = jl_contrasts(data), ...)
jlmer(
formula,
data,
family = NULL,
contrasts = jl_contrasts(data),
...,
fit = TRUE,
optsum = list(),
progress = interactive()
)Arguments
- formula
A formula written in Julia syntax. Can be a string or a language object.
- data
A data frame
- family
A distribution family
- contrasts
A Julia dictionary of contrasts Inferred from
databy default.- ...
Additional arguments to the
fit()function called in Julia- fit
Whether to fit the model. If
FALSE, returns the unfit model object.- optsum
A list of values to set for the optimizer. Special value
"debug"returns the OptSummary object for inspection.- progress
Whether to print model fitting progress. Defaults to
interactive()
Examples
# \donttest{
jlme_setup(restart = TRUE)
#> Starting Julia (v1.10.5) ...
#> Successfully set up Julia connection.
# Fixed effects models
lm(mpg ~ hp, mtcars)
#>
#> Call:
#> lm(formula = mpg ~ hp, data = mtcars)
#>
#> Coefficients:
#> (Intercept) hp
#> 30.09886 -0.06823
#>
jlm(mpg ~ hp, mtcars)
#> <Julia object of type StatsModels.TableRegressionModel>
#>
#> mpg ~ 1 + hp
#>
#> ────────────────────────────────────────────────────────────────────────────
#> Coef. Std. Error z Pr(>|z|) Lower 95% Upper 95%
#> ────────────────────────────────────────────────────────────────────────────
#> (Intercept) 30.0989 1.63392 18.42 <1e-75 26.8964 33.3013
#> hp -0.0682283 0.0101193 -6.74 <1e-10 -0.0880617 -0.0483948
#> ────────────────────────────────────────────────────────────────────────────
# Auto-handling of contrasts
x <- mtcars
x$cyl_helm <- factor(x$cyl)
contrasts(x$cyl_helm) <- contr.helmert(3)
colnames(contrasts(x$cyl_helm)) <- c("4vs6", "4&6vs8")
lm(mpg ~ cyl_helm, x)
#>
#> Call:
#> lm(formula = mpg ~ cyl_helm, data = x)
#>
#> Coefficients:
#> (Intercept) cyl_helm4vs6 cyl_helm4&6vs8
#> 20.502 -3.460 -2.701
#>
jlm(mpg ~ cyl_helm, x)
#> <Julia object of type StatsModels.TableRegressionModel>
#>
#> mpg ~ 1 + cyl_helm
#>
#> ─────────────────────────────────────────────────────────────────────────────
#> Coef. Std. Error z Pr(>|z|) Lower 95% Upper 95%
#> ─────────────────────────────────────────────────────────────────────────────
#> (Intercept) 20.5022 0.593528 34.54 <1e-99 19.3389 21.6655
#> cyl_helm: 4vs6 -3.46039 0.779174 -4.44 <1e-05 -4.98754 -1.93324
#> cyl_helm: 4&6vs8 -2.70108 0.387175 -6.98 <1e-11 -3.45993 -1.94223
#> ─────────────────────────────────────────────────────────────────────────────
# Mixed effects models
library(lme4)
glmer(r2 ~ Anger + Gender + (1 | id), VerbAgg, family = "binomial")
#> Generalized linear mixed model fit by maximum likelihood (Laplace
#> Approximation) [glmerMod]
#> Family: binomial ( logit )
#> Formula: r2 ~ Anger + Gender + (1 | id)
#> Data: VerbAgg
#> AIC BIC logLik deviance df.resid
#> 9504.505 9532.240 -4748.253 9496.505 7580
#> Random effects:
#> Groups Name Std.Dev.
#> id (Intercept) 1.059
#> Number of obs: 7584, groups: id, 316
#> Fixed Effects:
#> (Intercept) Anger GenderM
#> -1.10090 0.04626 0.26002
jlmer(r2 ~ Anger + Gender + (1 | id), VerbAgg, family = "binomial")
#> Error : Evaluation in Julia failed.
#> Original Julia error message:
#> KeyError: key 0x0000024f5a690e50 not found
#> Stacktrace:
#> [1] getindex
#> @ .\dict.jl:498 [inlined]
#> [2] decrefcount!(communicator::Main.RConnector.CommunicatoR{Sockets.TCPSocket}, ref::UInt64)
#> @ Main.RConnector C:\Users\jchoe\AppData\Local\R\win-library\4.4\JuliaConnectoR\Julia\sharing.jl:134
#> [3] decrefcounts(communicator::Main.RConnector.CommunicatoR{Sockets.TCPSocket}, refs::Vector{UInt8})
#> @ Main.RConnector C:\Users\jchoe\AppData\Local\R\win-library\4.4\JuliaConnectoR\Julia\sharing.jl:147
#> [4] invokelatest(::Any, ::Any, ::Vararg{Any}; kwargs::@Kwargs{})
#> @ Base .\essentials.jl:892
#> [5] invokelatest(::Any, ::Any, ::Vararg{Any})
#> @ Base .\essentials.jl:889
#> [6] evaluate!
#> @ C:\Users\jchoe\AppData\Local\R\win-library\4.4\JuliaConnectoR\Julia\evaluating.jl:67 [inlined]
#> [7] evaluate_checked!(call::Main.RConnector.Call, communicator::Main.RConnector.CommunicatoR{Sockets.TCPSocket})
#> @ Main.RConnector C:\Users\jchoe\AppData\Local\R\win-library\4.4\JuliaConnectoR\Julia\evaluating.jl:55
#> [8] serve_repl(sock::Sockets.TCPSocket)
#> @ Main.RConnector C:\Users\jchoe\AppData\Local\R\win-library\4.4\JuliaConnectoR\Julia\communicating.jl:185
#> [9] (::Main.RConnector.var"#3#5")()
#> @ Main.RConnector C:\Users\jchoe\AppData\Local\R\win-library\4.4\JuliaConnectoR\Julia\communicating.jl:145
#> <Julia object of type GeneralizedLinearMixedModel>
#>
#> r2 ~ 1 + Anger + Gender + (1 | id)
#>
#> Variance components:
#> Column VarianceStd.Dev.
#> id (Intercept) 1.12074 1.05865
#> ────────────────────────────────────────────────────
#> Coef. Std. Error z Pr(>|z|)
#> ────────────────────────────────────────────────────
#> (Intercept) -1.10115 0.280681 -3.92 <1e-04
#> Anger 0.0462741 0.0134906 3.43 0.0006
#> Gender: M 0.260057 0.153847 1.69 0.0910
#> ────────────────────────────────────────────────────
# Set optimizer options via `optsum`
jlmer(
r2 ~ Anger + Gender + (1 | id), VerbAgg, family = "binomial",
optsum = list(
optimizer = jl(":LN_NELDERMEAD"),
maxfeval = 10L
)
)
#> <Julia object of type GeneralizedLinearMixedModel>
#>
#> r2 ~ 1 + Anger + Gender + (1 | id)
#>
#> Variance components:
#> Column VarianceStd.Dev.
#> id (Intercept) 1.0 1.0
#> ────────────────────────────────────────────────────
#> Coef. Std. Error z Pr(>|z|)
#> ────────────────────────────────────────────────────
#> (Intercept) -0.851285 0.267259 -3.19 0.0014
#> Anger 0.0353915 0.012845 2.76 0.0059
#> Gender: M 0.204939 0.14652 1.40 0.1619
#> ────────────────────────────────────────────────────
stop_julia()
# }