R/predict.design.matrix.DDPstar.aux.R

In DDPstar: Density Regression via Dirichlet Process Mixtures of Normal Structured Additive Regression Models

predict_design.matrix.DDPstar.aux <-
function(object, newdata, select) {
	Xp <- NULL

	# Organize the newdataframe as it was in the original data
	cov.names <- names(object$iformula$data.cov)
	newdata <- newdata[, cov.names, drop = FALSE]

	# Standardised the continuous covariates
	newdata.std <- newdata
	cov.names.std <- colnames(object$iformula$cov.std)
	if(!is.null(cov.names.std)) {
		for(i in 1:length(cov.names.std)) {
			aux <- object$iformula$cov.std[,cov.names.std[i]]
			newdata.std[, cov.names.std[i]] <- (newdata[,cov.names.std[i]] - aux[1])/aux[2]
		}
	}

	for(i in select) {
		if(any(object$iformula$II[,i] == -1)) {
			if(object$iformula$h[i] == 0 | object$iformula$h[i] == 1) { # Linear and factor
				if(object$standardise) {
					mfp <- model.frame(object$terms[[i]], newdata.std, xlev = attr(object$terms[[i]], "xlev"))
				} else {
					mfp <- model.frame(object$terms[[i]], newdata, xlev = attr(object$terms[[i]], "xlev"))
				}
				Xp_aux <- model.matrix(object$terms[[i]], data = mfp, contrasts.arg = attr(object$terms[[i]], "contrast"))[,-1,drop = TRUE]
				Xp <- cbind(Xp, Xp_aux)
			} else if(object$iformula$h[i] == -1) { # Smooth effects
				Bs <- suppressWarnings(predict_bbase.bs(object$terms[[i]], newdata[,object$iformula$II[2,i], drop = TRUE]))
				Xp <- cbind(Xp, Bs)
			} else { # Random effects
				mfp <- model.frame(object$terms[[i]], newdata, xlev = attr(object$terms[[i]], "xlev"), na.action = na.pass)
				Xp_aux <- model.matrix(object$terms[[i]], data = mfp, contrasts.arg = attr(object$terms[[i]], "contrast"))
				Xp_aux <- Xp_aux[,-1,drop = FALSE]
				Xp_aux[is.na(Xp_aux)] <- 0
				Xp <- cbind(Xp, Xp_aux)
			}
		} else { # Factor by curve, varying coefficient or 2D
			if(object$iformula$by.var[i]) {
				if(is.factor(object$iformula$data.cov[,object$iformula$II[1,i]])) { # Factor by curve
					Bs <- predict_bbase.interaction.factor.by.curve.bs(object$terms[[i]], newdata[,object$iformula$II[2,i], drop = TRUE], newdata[,object$iformula$II[1,i], drop = TRUE])
					Xp <- cbind(Xp, Bs)
				} else { # Varying coefficient
					if(object$standardise) {
						Bs <- predict_bbase.interaction.vc.bs(object$terms[[i]], newdata[,object$iformula$II[2,i], drop = TRUE], newdata.std[,object$iformula$II[1,i], drop = TRUE])
					} else {
						Bs <- predict_bbase.interaction.vc.bs(object$terms[[i]], newdata[,object$iformula$II[2,i], drop = TRUE], newdata[,object$iformula$II[1,i], drop = TRUE])

					}
					Xp <- cbind(Xp, Bs)
				}
			} else { # 2D
				Bs <- predict_bbase.psanova.bs(object$terms[[i]], newdata[,object$iformula$II[1,i], drop = TRUE], newdata[,object$iformula$II[2,i], drop = TRUE])
				Xp <- cbind(Xp, Bs)
			}
		}
	}
	Xp
}