fit_blis: Fit Baseline-category Logit Intercept-Slope (BLIS) model on...
In ShinyItemAnalysis: Test and Item Analysis via Shiny
fit_blis R Documentation
Fit Baseline-category Logit Intercept-Slope (BLIS) model on nominal data
Description
blis fits the IRT Nominal Response Model to data from multiple-choice tests,
while accounting for the correct answer and treating this option as a baseline
in this baseline-category logit model. The intercept-slope parametrization in
BLIS can be converted to IRT (difficulty-discrimination) parametrization (BLIRT).
Usage
fit_blis(Data, key, ...)
blis(Data, key, ...)
Arguments
Data
data.frame or tibble with all columns being factors. Support
for matrix is limited and behavior not guaranteed.
key
A single-column data.frame, (not matrix) tibble or -
preferably - a factor vector of levels considered as correct responses.
...
Arguments passed on to mirt::mirt
SElogical; estimate the standard errors by computing the parameter information matrix?
See SE.type for the type of estimates available
covdataa data.frame of data used for latent regression models
formulaan R formula (or list of formulas) indicating how the latent traits
can be regressed using external covariates in covdata. If a named list
of formulas is supplied (where the names correspond to the latent trait names in model)
then specific regression effects can be estimated for each factor. Supplying a single formula
will estimate the regression parameters for all latent traits by default
itemdesigna data.frame with rows equal to the number of items and columns
containing any item-design effects. If items should be included in the design structure
(i.e., should be left in their canonical structure) then fewer rows can be used,
however the rownames must be defined and matched with colnames in the data
input. The item design matrix is constructed with the use of
item.formula. Providing this input will fix the associated 'd' intercepts
to 0, where applicable
item.formulaan R formula used to specify any intercept decomposition (e.g.,
the LLTM; Fischer, 1983). Note that only the right-hand side of the formula is required
for compensatory models.
For non-compensatory itemtypes (e.g., 'PC1PL') the formula must include
the name of the latent trait in the left hand side of the expression to indicate which
of the trait specification should have their intercepts decomposed (see MLTM; Embretson, 1984)
SE.typetype of estimation method to use for calculating the parameter information matrix
for computing standard errors and wald tests. Can be:
-
'Richardson', 'forward', or 'central' for the numerical Richardson,
forward difference, and central difference evaluation of observed Hessian matrix
-
'crossprod' and 'Louis' for standard error computations based on the variance of the
Fisher scores as well as Louis' (1982) exact computation of the observed information matrix.
Note that Louis' estimates can take a long time to obtain for large sample sizes and long tests
-
'sandwich' for the sandwich covariance estimate based on the
'crossprod' and 'Oakes' estimates (see Chalmers, 2018, for details)
-
'sandwich.Louis' for the sandwich covariance estimate based on the
'crossprod' and 'Louis' estimates
-
'Oakes' for Oakes' (1999) method using a central difference approximation
(see Chalmers, 2018, for details)
-
'SEM' for the supplemented EM (disables the accelerate option automatically; EM only)
-
'Fisher' for the expected information, 'complete' for information based
on the complete-data Hessian used in EM algorithm
-
'MHRM' and 'FMHRM' for stochastic approximations of observed information matrix
based on the Robbins-Monro filter or a fixed number of MHRM draws without the RM filter.
These are the only options supported when method = 'MHRM'
-
'numerical' to obtain the numerical estimate from a call to optim
when method = 'BL'
Note that both the 'SEM' method becomes very sensitive if the ML solution has
has not been reached with sufficient precision, and may be further sensitive
if the history of the EM cycles is not stable/sufficient for convergence of the respective estimates.
Increasing the number of iterations (increasing NCYCLES and decreasing
TOL, see below) will help to improve the accuracy, and can be
run in parallel if a mirtCluster object has been defined (this will be
used for Oakes' method as well). Additionally,
inspecting the symmetry of the ACOV matrix for convergence issues by passing
technical = list(symmetric = FALSE) can be helpful to determine if a sufficient
solution has been reached
methoda character object specifying the estimation algorithm to be used. The default is
'EM', for the standard EM algorithm with fixed quadrature, 'QMCEM' for
quasi-Monte Carlo EM estimation, or 'MCEM' for Monte Carlo EM estimation.
The option 'MHRM' may also be passed to use the MH-RM algorithm,
'SEM' for the Stochastic EM algorithm (first
two stages of the MH-RM stage using an optimizer other than a single Newton-Raphson iteration),
and 'BL' for the Bock and Lieberman
approach (generally not recommended for longer tests).
The 'EM' is generally effective with 1-3 factors, but methods such as the 'QMCEM',
'MCEM', 'SEM', or 'MHRM' should be used when the dimensions are 3 or more. Note that
when the optimizer is stochastic the associated SE.type is automatically changed to
SE.type = 'MHRM' by default to avoid the use of quadrature
optimizera character indicating which numerical optimizer to use. By default, the EM
algorithm will use the 'BFGS' when there are no upper and lower bounds box-constraints and
'nlminb' when there are.
Other options include the Newton-Raphson ('NR'),
which can be more efficient than the 'BFGS' but not as stable for more complex
IRT models (such as the nominal or nested logit models)
and the related 'NR1' which is also the Newton-Raphson
but consists of only 1 update that has been coupled with RM Hessian (only
applicable when the MH-RM algorithm is used). The MH-RM algorithm uses the 'NR1' by default,
though currently the 'BFGS', 'L-BFGS-B', and 'NR'
are also supported with this method (with
fewer iterations by default) to emulate stochastic EM updates.
As well, the 'Nelder-Mead' and 'SANN'
estimators are available, but their routine use generally is not required or recommended.
Additionally, estimation subroutines from the Rsolnp and nloptr
packages are available by passing the arguments 'solnp' and 'nloptr',
respectively. This should be used in conjunction with the solnp_args and
nloptr_args specified below. If equality constraints were specified in the
model definition only the parameter with the lowest parnum
in the pars = 'values' data.frame is used in the estimation vector passed
to the objective function, and group hyper-parameters are omitted.
Equality an inequality functions should be of the form function(p, optim_args),
where optim_args is a list of internally parameters that largely can be ignored
when defining constraints (though use of browser() here may be helpful)
dentypetype of density form to use for the latent trait parameters. Current options include
-
'Gaussian' (default) assumes a multivariate Gaussian distribution with an associated
mean vector and variance-covariance matrix
-
'empiricalhist' or 'EH' estimates latent distribution using an empirical histogram described by
Bock and Aitkin (1981). Only applicable for unidimensional models estimated with the EM algorithm.
For this option, the number of cycles, TOL, and quadpts are adjusted accommodate for
less precision during estimation (namely: TOL = 3e-5, NCYCLES = 2000, quadpts = 121)
-
'empiricalhist_Woods' or 'EHW' estimates latent distribution using an empirical histogram described by
Bock and Aitkin (1981), with the same specifications as in dentype = 'empiricalhist',
but with the extrapolation-interpolation method described by Woods (2007). NOTE: to improve stability
in the presence of extreme response styles (i.e., all highest or lowest in each item) the technical option
zeroExtreme = TRUE may be required to down-weight the contribution of these problematic patterns
-
'Davidian-#' estimates semi-parametric Davidian curves described by Woods and Lin (2009),
where the # placeholder represents the number of Davidian parameters to estimate
(e.g., 'Davidian-6' will estimate 6 smoothing parameters). By default, the number of
quadpts is increased to 121, and this method is only applicable for
unidimensional models estimated with the EM algorithm
Note that when itemtype = 'ULL' then a log-normal(0,1) density is used to support the unipolar scaling
constraina list of user declared equality constraints. To see how to define the
parameters correctly use pars = 'values' initially to see how the parameters are
labeled. To constrain parameters to be equal create a list with separate concatenated
vectors signifying which parameters to constrain. For example, to set parameters 1 and 5
equal, and also set parameters 2, 6, and 10 equal use
constrain = list(c(1,5), c(2,6,10)). Constraints can also be specified using the
mirt.model syntax (recommended)
calcNulllogical; calculate the Null model for additional fit statistics (e.g., TLI)?
Only applicable if the data contains no NA's and the data is not overly sparse
drawsthe number of Monte Carlo draws to estimate the log-likelihood for the MH-RM
algorithm. Default is 5000
survey.weightsa optional numeric vector of survey weights to apply for each case in the
data (EM estimation only). If not specified, all cases are weighted equally (the standard IRT
approach). The sum of the survey.weights must equal the total sample size for proper
weighting to be applied
quadptsnumber of quadrature points per dimension (must be larger than 2).
By default the number of quadrature uses the following scheme:
switch(as.character(nfact), '1'=61, '2'=31, '3'=15, '4'=9, '5'=7, 3).
However, if the method input is set to 'QMCEM' and this argument is left blank then
the default number of quasi-Monte Carlo integration nodes will be set to 5000 in total
TOLconvergence threshold for EM or MH-RM; defaults are .0001 and .001. If
SE.type = 'SEM' and this value is not specified, the default is set to 1e-5.
To evaluate the model using only the starting values pass TOL = NaN, and
to evaluate the starting values without the log-likelihood pass TOL = NA
gpcm_matsa list of matrices specifying how the scoring coefficients in the (generalized)
partial credit model should be constructed. If omitted, the standard gpcm format will be used
(i.e., seq(0, k, by = 1) for each trait). This input should be used if traits
should be scored different for each category (e.g., matrix(c(0:3, 1,0,0,0), 4, 2) for a
two-dimensional model where the first trait is scored like a gpcm, but the second trait is only
positively indicated when the first category is selected). Can be used when itemtypes
are 'gpcm' or 'Rasch', but only when the respective element in
gpcm_mats is not NULL
grsm.blockan optional numeric vector indicating where the blocking should occur when
using the grsm, NA represents items that do not belong to the grsm block (other items that may
be estimated in the test data). For example, to specify two blocks of 3 with a 2PL item for
the last item: grsm.block = c(rep(1,3), rep(2,3), NA). If NULL the all items are assumed
to be within the same group and therefore have the same number of item categories
rsm.blocksame as grsm.block, but for 'rsm' blocks
monopoly.ka vector of values (or a single value to repeated for each item) which indicate
the degree of the monotone polynomial fitted, where the monotone polynomial
corresponds to monopoly.k * 2 + 1 (e.g., monopoly.k = 2 fits a
5th degree polynomial). Default is monopoly.k = 1, which fits a 3rd degree polynomial
largea logical indicating whether unique response patterns should be obtained prior
to performing the estimation so as to avoid repeating computations on identical patterns.
The default TRUE provides the correct degrees of freedom for the model since all unique patterns
are tallied (typically only affects goodness of fit statistics such as G2, but also will influence
nested model comparison methods such as anova(mod1, mod2)), while FALSE will use the
number of rows in data as a placeholder for the total degrees of freedom. As such, model
objects should only be compared if all flags were set to TRUE or all were set to FALSE
Alternatively, if the collapse table of frequencies is desired for the purpose of saving computations
(i.e., only computing the collapsed frequencies for the data onte-time) then a character vector can
be passed with the arguement large = 'return' to return a list of all the desired
table information used by mirt. This list object can then be reused by passing it back
into the large argument to avoid re-tallying the data again
(again, useful when the dataset are very large and computing the tabulated data is
computationally burdensome). This strategy is shown below:
- Compute organized data
e.g., internaldat <- mirt(Science, 1, large = 'return')
- Pass the organized data to all estimation functions
e.g.,
mod <- mirt(Science, 1, large = internaldat)
GenRandomParslogical; generate random starting values prior to optimization instead of
using the fixed internal starting values?
acceleratea character vector indicating the type of acceleration to use. Default
is 'Ramsay', but may also be 'squarem' for the SQUAREM procedure (specifically,
the gSqS3 approach) described in Varadhan and Roldand (2008).
To disable the acceleration, pass 'none'
verboselogical; print observed- (EM) or complete-data (MHRM) log-likelihood
after each iteration cycle? Default is TRUE
solnp_argsa list of arguments to be passed to the solnp::solnp() function for
equality constraints, inequality constraints, etc
nloptr_argsa list of arguments to be passed to the nloptr::nloptr()
function for equality constraints, inequality constraints, etc
spline_argsa named list of lists containing information to be passed to the bs (default)
and ns for each spline itemtype. Each element must refer to the name of the itemtype with the
spline, while the internal list names refer to the arguments which are passed. For example, if item 2 were called
'read2', and item 5 were called 'read5', both of which were of itemtype 'spline' but item 5 should use the
ns form, then a modified list for each input might be of the form:
spline_args = list(read2 = list(degree = 4),
read5 = list(fun = 'ns', knots = c(-2, 2)))
This code input changes the bs() splines function to have a degree = 4 input,
while the second element changes to the ns() function with knots set a c(-2, 2)
controla list passed to the respective optimizers (i.e., optim(), nlminb(),
etc). Additional arguments have been included for the 'NR' optimizer: 'tol'
for the convergence tolerance in the M-step (default is TOL/1000), while the default
number of iterations for the Newton-Raphson optimizer is 50 (modified with the 'maxit'
control input)
technicala list containing lower level technical parameters for estimation. May be:
- NCYCLES
maximum number of EM or MH-RM cycles; defaults are 500 and 2000
- MAXQUAD
maximum number of quadratures, which you can increase if you have more than
4GB or RAM on your PC; default 20000
- theta_lim
range of integration grid for each dimension; default is c(-6, 6). Note that
when itemtype = 'ULL' a log-normal distribution is used and the range is change to
c(.01, and 6^2), where the second term is the square of the theta_lim input instead
- set.seed
seed number used during estimation. Default is 12345
- SEtol
standard error tolerance criteria for the S-EM and MHRM computation of the
information matrix. Default is 1e-3
- symmetric
logical; force S-EM/Oakes information matrix estimates to be symmetric? Default is TRUE
so that computation of standard errors are more stable. Setting this to FALSE can help
to detect solutions that have not reached the ML estimate
- SEM_window
ratio of values used to define the S-EM window based on the
observed likelihood differences across EM iterations. The default is
c(0, 1 - SEtol), which provides nearly the very full S-EM window (i.e.,
nearly all EM cycles used). To use the a smaller SEM window change the window to
to something like c(.9, .999) to start at a point farther into the EM history
- warn
logical; include warning messages during estimation? Default is TRUE
- message
logical; include general messages during estimation? Default is TRUE
- customK
a numeric vector used to explicitly declare the number of response
categories for each item. This should only be used when constructing mirt model for
reasons other than parameter estimation (such as to obtain factor scores), and requires
that the input data all have 0 as the lowest category. The format is the same as the
extract.mirt(mod, 'K') slot in all converged models
- customPriorFun
a custom function used to determine the normalized density for
integration in the EM algorithm. Must be of the form function(Theta, Etable){...},
and return a numeric vector with the same length as number of rows in Theta. The
Etable input contains the aggregated table generated from the current E-step
computations. For proper integration, the returned vector should sum to
1 (i.e., normalized). Note that if using the Etable it will be NULL
on the first call, therefore the prior will have to deal with this issue accordingly
- zeroExtreme
logical; assign extreme response patterns a survey.weight of 0
(formally equivalent to removing these data vectors during estimation)?
When dentype = 'EHW', where Woods' extrapolation is utilized,
this option may be required if the extrapolation causes expected densities to tend towards
positive or negative infinity. The default is FALSE
- customTheta
a custom Theta grid, in matrix form, used for integration.
If not defined, the grid is determined internally based on the number of quadpts
- nconstrain
same specification as the constrain list argument,
however imposes a negative equality constraint instead (e.g., a12 = -a21, which
is specified as nconstrain = list(c(12, 21))). Note that each specification
in the list must be of length 2, where the second element is taken to be -1 times the
first element
- delta
the deviation term used in numerical estimates when computing the ACOV matrix
with the 'forward' or 'central' numerical approaches, as well as Oakes' method with the
Richardson extrapolation. Default is 1e-5
- parallel
logical; use the parallel cluster defined by mirtCluster?
Default is TRUE
- storeEMhistory
logical; store the iteration history when using the EM algorithm?
Default is FALSE. When TRUE, use extract.mirt to extract
- internal_constraints
logical; include the internal constraints when using certain
IRT models (e.g., 'grsm' itemtype). Disable this if you want to use special optimizers
such as the solnp. Default is TRUE
- gain
a vector of two values specifying the numerator and exponent
values for the RM gain function (val1 / cycle)^val2.
Default is c(0.10, 0.75)
- BURNIN
number of burn in cycles (stage 1) in MH-RM; default is 150
- SEMCYCLES
number of SEM cycles (stage 2) in MH-RM; default is 100
- MHDRAWS
number of Metropolis-Hasting draws to use in the MH-RM at each iteration; default is 5
- MHcand
a vector of values used to tune the MH sampler. Larger values will
cause the acceptance ratio to decrease. One value is required for each group in
unconditional item factor analysis (mixedmirt() requires additional values
for random effect). If null, these values are determined internally, attempting to
tune the acceptance of the draws to be between .1 and .4
- MHRM_SE_draws
number of fixed draws to use when SE=TRUE and SE.type = 'FMHRM'
and the maximum number of draws when SE.type = 'MHRM'. Default is 2000
- MCEM_draws
a function used to determine the number of quadrature points to draw for the
'MCEM' method. Must include one argument which indicates the iteration number of the
EM cycle. Default is function(cycles) 500 + (cycles - 1)*2, which starts the number of
draws at 500 and increases by 2 after each full EM iteration
- info_if_converged
logical; compute the information matrix when using the MH-RM algorithm
only if the model converged within a suitable number of iterations? Default is TRUE
- logLik_if_converged
logical; compute the observed log-likelihood when using the MH-RM algorithm
only if the model converged within a suitable number of iterations? Default is TRUE
- keep_vcov_PD
logical; attempt to keep the variance-covariance matrix of the latent traits
positive definite during estimation in the EM algorithm? This generally improves the convergence
properties when the traits are highly correlated. Default is TRUE
Details
For the details on coef method dispatched for fitted BLIS model, see
coef,BlisClass-method. To get more on the class, see BlisClass.
Value
Fitted model of class BlisClass (extending standard mirt's
SingleGroupClass).
Author(s)
Jan Netik
Institute of Computer Science of the Czech Academy of
Sciences
netik@cs.cas.cz
Patricia Martinkova
Institute of Computer Science of the Czech Academy
of Sciences
martinkova@cs.cas.cz
See Also
Other BLIS/BLIRT related:
BlisClass-class,
coef,BlisClass-method,
get_orig_levels(),
nominal_to_int(),
obtain_nrm_def(),
print.blis_coefs()
Examples
fitted_blis <- fit_blis(HCItest[, 1:20], HCIkey, SE = TRUE)
coef(fitted_blis)
coef(fitted_blis)$`Item 12`
coef(fitted_blis, IRTpars = TRUE)
coef(fitted_blis, IRTpars = TRUE, CI = 0.90) # 90% CI instead of 95% CI
coef(fitted_blis, IRTpars = TRUE, printSE = TRUE) # SE instead of CI
ShinyItemAnalysis documentation built on Aug. 8, 2025, 6:40 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
| fit_blis | R Documentation |
Fit Baseline-category Logit Intercept-Slope (BLIS) model on nominal data
Description
blis fits the IRT Nominal Response Model to data from multiple-choice tests,
while accounting for the correct answer and treating this option as a baseline
in this baseline-category logit model. The intercept-slope parametrization in
BLIS can be converted to IRT (difficulty-discrimination) parametrization (BLIRT).
Usage
fit_blis(Data, key, ...)
blis(Data, key, ...)
Arguments
Data |
data.frame or tibble with all columns being factors. Support for matrix is limited and behavior not guaranteed. |
key |
A single-column |
... |
Arguments passed on to
|
Details
For the details on coef method dispatched for fitted BLIS model, see
coef,BlisClass-method. To get more on the class, see BlisClass.
Value
Fitted model of class BlisClass (extending standard mirt's
SingleGroupClass).
Author(s)
Jan Netik
Institute of Computer Science of the Czech Academy of
Sciences
netik@cs.cas.cz
Patricia Martinkova
Institute of Computer Science of the Czech Academy
of Sciences
martinkova@cs.cas.cz
See Also
Other BLIS/BLIRT related:
BlisClass-class,
coef,BlisClass-method,
get_orig_levels(),
nominal_to_int(),
obtain_nrm_def(),
print.blis_coefs()
Examples
fitted_blis <- fit_blis(HCItest[, 1:20], HCIkey, SE = TRUE)
coef(fitted_blis)
coef(fitted_blis)$`Item 12`
coef(fitted_blis, IRTpars = TRUE)
coef(fitted_blis, IRTpars = TRUE, CI = 0.90) # 90% CI instead of 95% CI
coef(fitted_blis, IRTpars = TRUE, printSE = TRUE) # SE instead of CI
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.