lavaan: Latent Variable Analysis

# export to Mplus syntax

lav2mplus <- function(lav, group.label = NULL) {
  lav <- lav2check(lav)
  header <- "  ! this model syntax is autogenerated by lavExport\n"
  footer <- "\n"

  lav <- as.data.frame(lav, stringsAsFactors = FALSE)
  ngroups <- lav_partable_ngroups(lav)

  lav_one_group <- function(lav) {
    # mplus does not like variable names with a 'dot'
    # replace them by an underscore '_'
    lav$lhs <- gsub("\\.", "_", lav$lhs)
    lav$rhs <- gsub("\\.", "_", lav$rhs)

    # remove contraints (:=, <, >, ==) here
    con.idx <- which(lav$op %in% c(":=", "<", ">", "=="))
    if (length(con.idx) > 0L) {
      lav <- lav[-con.idx, ]
    }

    # remove exogenous variances/covariances/intercepts...
    exo.idx <- which(lav$exo == 1L & lav$op %in% c("~~", "~1"))
    if (length(exo.idx)) {
      lav <- lav[-exo.idx, ]
    }

    # remove intercepts for categorical variables
    ord.names <- unique(lav$lhs[lav$op == "|"])
    ord.int.idx <- which(lav$op == "~1" & lav$lhs %in% ord.names)
    if (length(ord.int.idx)) {
      lav <- lav[-ord.int.idx, ]
    }

    # end of line
    lav$eol <- rep(";", length(lav$lhs))
    lav$ustart <- ifelse(is.na(lav$ustart), "", lav$ustart)
    lav$rhs2 <- ifelse(lav$free == 0L,
      paste("@", lav$ustart, sep = ""),
      paste("*", lav$ustart, sep = "")
    )
    lav$plabel <- gsub("\\.", "", lav$plabel)
    LABEL <- ifelse(lav$label == "", lav$plabel, lav$label)
    lav$plabel <- ifelse(LABEL == "", LABEL,
      paste(" (", LABEL, ")", sep = "")
    )

    # remove variances for ordered variables
    ov.names.ord <- vnames(lav, type = "ov.ord")
    ord.idx <- which(lav$lhs %in% ov.names.ord &
      lav$op == "~~" &
      lav$free == 0L &
      lav$lhs == lav$rhs)
    lav$lhs[ord.idx] <- paste("! ", lav$lhs[ord.idx], sep = "")
    lav$op[ord.idx] <- ""
    lav$rhs[ord.idx] <- ""

    # variances
    var.idx <- which(lav$op == "~~" & lav$rhs == lav$lhs)
    lav$op[var.idx] <- ""
    lav$rhs[var.idx] <- ""

    # scaling factors
    scal.idx <- which(lav$op == "~*~")
    lav$op[scal.idx] <- ""
    lav$rhs2[scal.idx] <- paste(lav$rhs2[scal.idx], "}", sep = "")
    lav$lhs[scal.idx] <- "{"

    # intercepts - excluding categorical observed
    int.idx <- which(lav$op == "~1")
    lav$op[int.idx] <- ""
    lav$rhs2[int.idx] <- paste(lav$rhs2[int.idx], "]", sep = "")
    lav$lhs[int.idx] <- paste("[", lav$lhs[int.idx], sep = "")

    # thresholds
    th.idx <- which(lav$op == "|")
    lav$op[th.idx] <- "$"
    lav$rhs[th.idx] <- gsub("t", "", x = lav$rhs[th.idx])
    lav$rhs2[th.idx] <- paste(lav$rhs2[th.idx], "]", sep = "")
    lav$lhs[th.idx] <- paste("[", lav$lhs[th.idx], sep = "")

    # replace binary operators
    lav$op <- ifelse(lav$op == "=~", " BY ", lav$op)
    lav$op <- ifelse(lav$op == "~", " ON ", lav$op)
    lav$op <- ifelse(lav$op == "~~", " WITH ", lav$op)


    lav2 <- paste(lav$lhs, lav$op, lav$rhs, lav$rhs2,
      lav$plabel, lav$eol,
      sep = ""
    )

    body <- paste(" ", lav2, collapse = "\n")

    body
  }

  if (ngroups == 1L) {
    body <- lav_one_group(lav)
  } else {
    group.values <- lav_partable_group_values(lav)
    # group 1
    body <- lav_one_group(lav[lav$group == group.values[1], ])

    if (is.null(group.label) || length(group.label) == 0L) {
      group.label <- paste(1:ngroups)
    }

    for (g in 2:ngroups) {
      body <- paste(body,
        paste("\nMODEL ", group.label[g], ":\n", sep = ""),
        lav_one_group(lav[lav$group == group.values[g], ]),
        sep = ""
      )
    }
  }

  # constraints go to a 'MODEL CONSTRAINTS' block
  con.idx <- which(lav$op %in% c(":=", "<", ">", "=="))
  if (length(con.idx) > 0L) {
    ### FIXME: we need to convert the operator
    ### eg b^2 --> b**2, others??
    lav$lhs[con.idx] <- gsub("\\^", "**", lav$lhs[con.idx])
    lav$rhs[con.idx] <- gsub("\\^", "**", lav$rhs[con.idx])

    constraints <- "\nMODEL CONSTRAINT:\n"
    # define 'new' variables
    def.idx <- which(lav$op == ":=")
    if (length(def.idx) > 0L) {
      def <- paste(lav$lhs[def.idx], collapse = " ")
      constraints <- paste(constraints, "NEW (", def, ");")
      lav$op[def.idx] <- "="
    }
    # replace '==' by '='
    eq.idx <- which(lav$op == "==")
    if (length(eq.idx) > 0L) {
      lav$op[eq.idx] <- "="
    }
    con <- paste(gsub("\\.", "", lav$lhs[con.idx]), " ",
      lav$op[con.idx], " ",
      gsub("\\.", "", lav$rhs[con.idx]), ";",
      sep = ""
    )
    con2 <- paste("  ", con, collapse = "\n")
    constraints <- paste(constraints, con2, sep = "\n")
  } else {
    constraints <- ""
  }

  out <- paste(header, body, constraints, footer, sep = "")
  class(out) <- c("lavaan.character", "character")
  out
}


# helper functions
lav_mplus_estimator <- function(object) {
  estimator <- object@Options$estimator
  if (estimator == "DWLS") {
    estimator <- "WLS"
  }

  # only 1 argument for 'test' is allowed
  if (length(object@Options$test) > 1L) {
    standard.idx <- which(object@Options$test == "standard")
    if (length(standard.idx) > 1L) {
      object@Options$test <- object@Options$test[-standard.idx]
    }
    if (length(object@Options$test) > 1L) {
      lav_msg_warn(gettext("only first (non-standard) test will be used"))
      object@Options$test <- object@Options$test[1]
    }
  }

  if (estimator == "ML") {
    if (object@Options$test %in% c("yuan.bentler", "yuan.bentler.mplus")) {
      estimator <- "MLR"
    } else if (object@Options$test == "satorra.bentler") {
      estimator <- "MLM"
    } else if (object@Options$test == "scaled.shifted") {
      estimator <- "MLMV"
    } else if (object@Options$se == "first.order") {
      estimator <- "MLF"
    }
  } else if (estimator %in% c("ULS", "WLS")) {
    if (object@Options$test == "satorra.bentler") {
      estimator <- paste(estimator, "M", sep = "")
    } else if (object@Options$test == "scaled.shifted") {
      estimator <- paste(estimator, "MV", sep = "")
    }
  } else if (estimator == "MML") {
    estimator <- "ML"
  }

  estimator
}

lav_mplus_header <- function(data.file = NULL, group.label = "", ov.names = "",
                             listwise = FALSE,
                             ov.ord.names = "", estimator = "ML",
                             meanstructure = FALSE,
                             weight.name = character(0L),
                             information = "observed",
                             data.type = "full", nobs = NULL) {
  # replace '.' by '_' in all variable names
  ov.names <- gsub("\\.", "_", ov.names)
  ov.ord.names <- gsub("\\.", "_", ov.ord.names)

  ### FIXME!!
  ### this is old code from lavaan 0.3-1
  ### surely, this can be done better...

  # TITLE command
  c.TITLE <- "TITLE:\n"
  c.TITLE <- paste(
    c.TITLE,
    "  [This syntax is autogenerated by lavExport]\n"
  )

  # DATA command
  c.DATA <- "DATA:\n"
  ngroups <- length(data.file)
  if (ngroups == 1L) {
    c.DATA <- paste(c.DATA,
      "  file is ", data.file, ";\n",
      sep = ""
    )
  } else {
    for (g in 1:ngroups) {
      c.DATA <- paste(c.DATA,
        "  file (", group.label[g], ") is ",
        data.file[g], ";\n",
        sep = ""
      )
    }
  }
  if (data.type == "full") {
    c.DATA <- paste(c.DATA, "  type is individual;\n", sep = "")
    if (listwise) {
      c.DATA <- paste(c.DATA, "  listwise = on;\n", sep = "")
    }
  } else if (data.type == "moment") {
    c.DATA <- paste(c.DATA, "  type is fullcov;\n", sep = "")
    c.DATA <- paste(c.DATA, "  nobservations are ", nobs, ";\n", sep = "")
  } else {
    lav_msg_stop(gettext("data.type must be full or moment"))
  }

  # VARIABLE command
  c.VARIABLE <- "VARIABLE:\n"
  c.VARIABLE <- paste(c.VARIABLE, "  names are", sep = "")
  nvar <- length(ov.names)
  tmp <- 0
  for (i in 1:nvar) {
    if (tmp %% 6 == 0) {
      c.VARIABLE <- paste(c.VARIABLE, "\n    ", sep = "")
    }
    c.VARIABLE <- paste(c.VARIABLE, ov.names[i], sep = " ")
    tmp <- tmp + 1
  }
  c.VARIABLE <- paste(c.VARIABLE, ";\n", sep = "")
  # missing
  if (data.type == "full") {
    c.VARIABLE <- paste(c.VARIABLE,
      "  missing are all (-999999);\n",
      sep = ""
    )
  }
  # categorical?
  if (length(ov.ord.names)) {
    c.VARIABLE <- paste(c.VARIABLE, "  categorical are", sep = "")
    nvar <- length(ov.ord.names)
    tmp <- 0
    for (i in 1:nvar) {
      if (tmp %% 6 == 0) {
        c.VARIABLE <- paste(c.VARIABLE, "\n    ", sep = "")
      }
      c.VARIABLE <- paste(c.VARIABLE, ov.ord.names[i])
      tmp <- tmp + 1
    }
    c.VARIABLE <- paste(c.VARIABLE, ";\n", sep = "")
  }
  # weight variable?
  if (length(weight.name) > 0L) {
    c.VARIABLE <- paste(c.VARIABLE,
      "  weight = ", weight.name, ";\n",
      sep = ""
    )
  }

  # ANALYSIS command
  c.ANALYSIS <- paste("ANALYSIS:\n  type = general;\n", sep = "")
  c.ANALYSIS <- paste(c.ANALYSIS, "  estimator = ", toupper(estimator),
    ";\n",
    sep = ""
  )
  if (toupper(estimator) %in% c("ML", "MLR")) {
    c.ANALYSIS <- paste(c.ANALYSIS, "  information = ", information[1],
      ";\n",
      sep = ""
    )
  }
  if (!meanstructure) {
    c.ANALYSIS <- paste(c.ANALYSIS, "  model = nomeanstructure;\n",
      sep = ""
    )
  }

  # MODEL command
  c.MODEL <- paste("MODEL:\n")

  # assemble pre-model header
  out <- paste(c.TITLE, c.DATA, c.VARIABLE, c.ANALYSIS, c.MODEL, sep = "")

  out
}