umxACE_DE: Build a Double-Entry Censored Data ACE twin model
In umx: Structural Equation Modeling and Twin Modeling in R

umxACE_DE

R Documentation

Build a Double-Entry Censored Data ACE twin model

Description

umxACE_DE implements a twin model for handling censored/floor-effect variables using a "double-entry" strategy. Each censored variable is represented by a pair of columns: one continuous (holding non-censored values) and one binary/ordinal (holding the censored indicator). The path coefficients (A, C, E) are constrained to be equal across the two paired columns.

Fully continuous variables may be mixed with double-entry pairs in selDVs (e.g. c("ht", "wt_cont", "wt_cens")). Prepare the censored traits with umx_make_double_entry_data(). At least one contiguous ⁠_cont⁠/⁠_cens⁠ pair is required; for all-continuous models use umxACE().

"Trust the Science" note

If you know the censoring value, keep fixCensorThresholds = "yes" / "auto" to lock that in. Estimated-free threshold can "wander" wildly. Reserve it only for unknown cut, and even then treat estimated thresholds with suspicion. Free threshold will "run" but it is not great at recovering the true floor/threshold.

Usage

umxACE_DE(
  name = "ACE_DE",
  selDVs,
  selCovs = NULL,
  dzData = NULL,
  mzData = NULL,
  sep = "_T",
  data = NULL,
  zyg = "zygosity",
  fixCensorThresholds = c("yes", "auto", "no"),
  censorCuts = NULL,
  doubleEntrySuffix = c("_cont", "_cens"),
  type = c("Auto", "FIML", "cov", "cor", "WLS", "DWLS", "ULS"),
  numObsDZ = NULL,
  numObsMZ = NULL,
  boundDiag = 0,
  allContinuousMethod = c("cumulants", "marginals"),
  autoRun = getOption("umx_auto_run"),
  intervals = FALSE,
  tryHard = c("no", "yes", "ordinal", "search"),
  optimizer = NULL,
  nSib = 2,
  dzAr = 0.5,
  dzCr = 1,
  weightVar = NULL,
  equateMeans = TRUE,
  addStd = TRUE,
  addCI = TRUE
)

Arguments

`name`	The name of the model (defaults to "ACE").
`selDVs`	Base names of variables to model. Include fully observed continuous traits by base name (e.g. `"ht"`). Include each censored trait as an adjacent pair of prepped names (e.g. `⁠"wt_cont", "wt_cens"⁠`). Prep censored data with `umx_make_double_entry_data()` first.
`selCovs`	(optional) covariates to include from the data (do not include sep in names)
`dzData`	The DZ dataframe.
`mzData`	The MZ dataframe.
`sep`	The separator in twin variable names, often "_T", e.g. "dep_T1".
`data`	If provided, dzData and mzData are treated as levels of zyg to select() MZ and DZ data sets (default = NULL)
`zyg`	If data provided, this column is used to select rows by zygosity (Default = "zygosity")
`fixCensorThresholds`	One of `"yes"` (default), `"auto"`, or `"no"`. `"yes"`: fix every double-entry pair from `censorCuts` and/or prep attribute `umxDoubleEntry`. `"auto"`: fix only pairs with a finite known cut in `censorCuts` or the prep attribute. `"no"`: free binary thresholds.
`censorCuts`	Optional named numeric vector of known cuts on the analysis scale (after any scaling). Names may be trait base (`"wt"`), continuous base (`"wt_cont"`), or censored base (`"wt_cens"`). When non-`NULL`, fixCensorThresholds must be `"yes"` or `"auto"`; only the named pairs are fixed (partial `censorCuts` is allowed).
`doubleEntrySuffix`	Suffixes for the continuous and censored variables (default = c("_cont", "_cens")).
`type`	Analysis method one of c("Auto", "FIML", "cov", "cor", "WLS", "DWLS", "ULS")
`numObsDZ`	Number of DZ twins: Set this if you input covariance data.
`numObsMZ`	Number of MZ twins: Set this if you input covariance data.
`boundDiag`	Numeric lbound for diagonal of the a, c, and e matrices. Defaults to 0.
`allContinuousMethod`	"cumulants" or "marginals". Used in all-continuous WLS data to determine if a means model needed.
`autoRun`	Whether to run the model (default), or just to create it and return without running.
`intervals`	Whether to run mxCI confidence intervals (default = FALSE)
`tryHard`	Default ('no') uses normal mxRun. "yes" uses mxTryHard. Other options: "ordinal", "search"
`optimizer`	Optionally set the optimizer (default NULL does nothing).
`nSib`	Number of siblings in a family (default = 2).
`dzAr`	The DZ genetic correlation (defaults to .5, vary to examine assortative mating).
`dzCr`	The DZ "C" correlation (defaults to 1: set to .25 to make an ADE model).
`weightVar`	If provided, a vector objective will be used to weight the data. (default = NULL).
`equateMeans`	Whether to equate the means across twins (defaults to TRUE).
`addStd`	Whether to add the algebras to compute a std model (defaults to TRUE).
`addCI`	Whether to add intervals to compute CIs (defaults to TRUE).

Details

Double-entry modeling represents a floor- or ceiling-censored trait using an adjacent pair of manifest columns: one continuous (⁠_cont⁠) and one binary/ordinal factor (⁠_cens⁠).

Likelihood Evaluation & Missingness Structure: Each individual contributes exactly one non-missing likelihood element for the censored trait:

Non-censored cases (x > cut): ⁠_cont⁠ contains the observed numeric value evaluated via the continuous normal density f(x), while ⁠_cens⁠ is set to NA.
Censored cases (x \le cut): ⁠_cont⁠ is set to NA, while ⁠_cens⁠ contains the ordinal factor level evaluated via the cumulative threshold probability P(Y \le \tau).

Leaving ⁠_cens⁠ non-missing for observed continuous rows would double-count the tail density (evaluating both f(x) and P(Y > \tau) for the same individual), introducing artificial covariance dependencies and inflating density estimates.

Value

OpenMx::mxModel() of subclass mxModel.ACE

Examples

require(umx)

##################################
# Toy example 1: Threshold known #
##################################

# Weighing scales lowest value is 60kg, with lower values registering as 60kg.

# NOTE: We first put height on a scale comparable to weight to ease estimation
data(twinData)
twinData[, c("ht1", "ht2")] = twinData[, c("ht1", "ht2")] *100 # metre->cm

clinic = twinData
clinic$wt1[!is.na(twinData$wt1) & (twinData$wt1 <= 60)] = 60
clinic$wt2[!is.na(twinData$wt2) & (twinData$wt2 <= 60)] = 60

# Double-entry prep for weight (floor at 0) Creates "wt_cont*" and "wt_cens*"
prep = umx_make_double_entry_data(clinic, cols = list(wt = "<= 60"), sep = "")

mzData = prep[prep$zygosity %in% "MZFF", ]
dzData = prep[prep$zygosity %in% "DZFF", ]


# 1. Correct mixed model: continuous height + double-entry censored weight
# Default: free thresholds. For known LOD, use fixCensorThresholds = "auto" or "yes"
# with censorCuts = c(wt = 0) (cut on analysis scale; means cont=cens equated).

mDE = umxACE_DE(name = "htWtDE", mzData = mzData, dzData = dzData, sep = "",
	selDVs = c("ht", "wt_cont", "wt_cens"), tryHard = "yes")

# Table: Model Fit Summary for 'htWtDE'
# 
# |Model  | EP|  -2LL|   df|   AIC|      BIC|
# |:------|--:|-----:|----:|-----:|--------:|
# |htWtDE | 11| 37082| 7781| 37104| 37165.51|


# Table: Standardized parameter estimates from a 2-trait double-entry Cholesky ACE model.
# A: additive genetic; C: common environment; E: unique environment.
# 
# |        |    a1|a2    |    c1|c2    |     e1|e2    |
# |:-------|-----:|:-----|-----:|:-----|------:|:-----|
# |ht      | 0.874|      | 0.373|      |  0.312|      |
# |wt_cens | 0.437|0.381 | 0.079|0.301 | -0.148|0.738 |
# Double-entry thresholds fixed: wt_cens @ 60 (means equated to wt_cont).


# Table: Means (from model$top$expMean)
# 
# |          |     ht1| wt_cont1| wt_cens1|     ht2| wt_cont2| wt_cens2|
# |:---------|-------:|--------:|--------:|-------:|--------:|--------:|
# |intercept | 162.499|   56.961|   56.961| 162.499|   56.961|   56.961|


####################################
# Toy example 2: Threshold unknown #
####################################

# GP only gets paid to record weight for higher-BMI people (~BMI 22 here).
# Others are coded "0").
# NOTE: We don't know what threshold 0 mapped to, so we *estimate* the threshold.
# NOTE: This example is contrived, and creates a nasty collider with height.

# Height fully observed; Weighing scales lowest value is 60kg, with all
# values beneath this registering as 60kg.

# NOTE: First we get the height data on a comparable scale to weight to ease model estimation
data(twinData)
twinData[, c("ht1", "ht2")] = twinData[, c("ht1", "ht2")] * 100 # metres -> cms

bmiCut = 22 # about the 80th percentile in this sample.

clinic = twinData
clinic$wt1[!is.na(twinData$wt1) & (twinData$bmi1 <= 22)] = 0
clinic$wt2[!is.na(twinData$wt2) & (twinData$bmi2 <= 22)] = 0

# Double-entry prep for weight (floor at 0) Creates "wt_cont*" and "wt_cens*"
prep = umx_make_double_entry_data(data = clinic, cols = list(wt = 0), sep = "")

mzData = prep[prep$zygosity %in% "MZFF", ]
dzData = prep[prep$zygosity %in% "DZFF", ]

# 1. Correct mixed model: continuous height + double-entry censored weight
# Default: free thresholds. For known LOD, use fixCensorThresholds = "auto" or "yes"
# with censorCuts = c(wt = 0) (cut on analysis scale; means cont=cens equated).

# Note how we set fixCensorThresholds = "no"

mDE = umxACE_DE(name = "htWtDE", mzData = mzData, dzData = dzData, sep = "",
	selDVs = c("ht", "wt_cont", "wt_cens"), fixCensorThresholds="no", tryHard = "yes")

# Table: Model Fit Summary for 'htWtDE'
# 
# |Model  | EP|     -2LL|   df|      AIC|      BIC|
# |:------|--:|--------:|----:|--------:|--------:|
# |htWtDE | 12| 32406.96| 7780| 32430.96| 32498.06|
# 
# 
# Table: Standardized parameter estimates from a 2-trait double-entry Cholesky ACE model. 
#   A: additive genetic; C: common environment; E: unique environment.
# 
# |        |    a1|a2    |    c1|c2    |    e1|e2    |
# |:-------|-----:|:-----|-----:|:-----|-----:|:-----|
# |ht      | 0.883|      | 0.268|      | 0.385|      |
# |wt_cens | 0.228|0.105 | 0.628|0.235 | 0.180|0.675 |
# 
# 
# Table: Means (from model$top$expMean)
# 
# |          |     ht1| wt_cont1|wt_cens1 |     ht2| wt_cont2|wt_cens2 |
# |:---------|-------:|--------:|:--------|-------:|--------:|:--------|
# |intercept | 162.918|   68.623|0        | 162.918|   68.623|0        |

umxSummary(mDE, std = TRUE)

# 2. Gold standard: uncensored bivariate ACE on true height and weight
mzTrue = twinData[twinData$zygosity %in% "MZFF", ]
dzTrue = twinData[twinData$zygosity %in% "DZFF", ]
mTrue  = umxACE("htWtTrue", selDVs = c("ht", "wt"), 
           mzData = mzTrue, dzData = dzTrue, sep = "", tryHard = "yes")

# Table: Standardized parameter estimates from 2-factor Cholesky ACE
#   A: additive genetic; C: common environment; E: unique environment.
# |   |    a1|a2    |     c1|c2 |    e1|e2    |
# |:--|-----:|:-----|------:|:--|-----:|:-----|
# |ht | 0.899|      |  0.252|   | 0.357|      |
# |wt | 0.419|0.766 | -0.049|.  | 0.074|0.479 |

umxSummary(mTrue, std = TRUE)

# 3. Mistaken analysis: treat 0 (not weighed) as a real continuous weight
mzNaive = clinic[clinic$zygosity %in% "MZFF", ]
dzNaive = clinic[clinic$zygosity %in% "DZFF", ]
mNaive = umxACE("htWtNaive0", selDVs= c("ht", "wt"), 
            mzData= mzNaive, dzData= dzNaive, sep ="",tryHard = "yes")
umxSummary(mNaive, std = TRUE)

# Naive mean(wt) pulled toward 0, covariance obscured; prefer mDE when zeros mean censored.
#
# Table: Standardized parameter estimates from a 2-factor Cholesky ACE model.
#   A: additive genetic; C: common environment; E: unique environment.
#
# |   |    a1|a2    |     c1|c2 |     e1|e2    |
# |:--|-----:|:-----|------:|:--|------:|:-----|
# |ht | 0.899|      |  0.252|   |  0.357|      |
# |wt | 0.057|0.717 | -0.238|.  | -0.072|0.649 |
#
#
# Table: Means (from model$top$expMean)
# 
# |          |     ht1|    wt1|     ht2|    wt2|
# |:---------|-------:|------:|-------:|------:|
# |intercept | 162.515| 20.984| 162.515| 20.984|

umx documentation built on Oct. 11, 2026, 5:07 p.m.