faster alternatives to sparse.model.matrix?
Try that but remove LHS variables from dataset before hands? We use that for the data processing part of putting our xgboost models in production.

fspmatrix &lt;- function(df) {
  # Building sparse matrix with sparseMatrix to save time vs sparse.model.matrix

  # Details about the df to help build triplet
  nr &lt;- nrow(df)
  nc &lt;- ncol(df)
  nlevels &lt;- pmax(1, sapply(df, nlevels))
  fac &lt;- which(sapply(df, is.factor))
  # Check that all factors have at least 2 levels
  stopifnot(min(nlevels[fac]) &gt; 1)

  # Building the value vector removing zero positions and first levels of factors after the first factor
  x &lt;- as.numeric(unlist(unlist(df)))
  jx &lt;- x # for use in building j
  # Positions of values coming from factors
  pos &lt;- as.vector(sapply(fac * nr, function(x) x - nr + (1:nr)))
  # Positions of first levels of factors after the first factor
  pos1 &lt;- pos[-1:-nr]
  pos1 &lt;- pos1[which(x[pos1] == 1)]
  # Replacing the values of factors with 1 for one hot encoding
  x[pos] &lt;- 1
  # Positions of non-zero values
  nzpos &lt;- which(x[-pos1] != 0)
  # Values to build triplet around
  x &lt;- x[-pos1][nzpos]

  # Building the i vector, the first dim part of the triplet
  i &lt;- rep(seq_len(nr), nc)[-pos1][nzpos]

  # Building the j vector, the second dim part of the triplet
  jnlevels &lt;- nlevels
  jnlevels[fac[-1]] &lt;- jnlevels[fac[-1]] - 1
  jx[-pos] &lt;- 1
  jx[pos[-1:-nr]] &lt;- jx[pos[-1:-nr]] - 1
  cs &lt;- cumsum(c(0, head(jnlevels, -1)))
  j &lt;- jx + rep(cs, each = nr)
  j &lt;- j[-pos1][nzpos]

  # Building the dimnames vector
  tag &lt;- TRUE
  fspnames &lt;- as.vector(unlist(sapply(names(df), function(n) {
    clvl &lt;- as.character(levels(df[[n]]))
    isL &lt;- is.logical(df[[n]])
    isF &lt;- is.factor(df[[n]])
    indx &lt;- if (tag &amp; isF) {
      1:length(clvl)
    } else {
      -1
    }
    if (tag &amp; isF) tag &lt;&lt;- FALSE
    paste0(n, clvl[indx], ifelse(isL, "TRUE", ""))
  })))

  fspm &lt;- sparseMatrix(i = i, j = j, x = x, dims = c(nr, tail(cs, 1) + tail(jnlevels, 1)), dimnames = list(NULL, fspnames))

  return(fspm)
}