SCM

[#2126] predict for censReg

View Trackers | Feature Requests | Download .csv | Monitor

Date:
2012-07-05 12:38
Priority:
3
State:
Open
Submitted by:
Tomas Krehlik (tomaskrehlik)
Assigned to:
Arne Henningsen (arne)
Product:
censReg
Operating System:
All
Status:
 
Summary:
predict for censReg

Detailed description
Hello,

I pretty miss some predict functionality in censReg, that would at provide fitted values. And to be constructive, I am appending code, which seems to be pretty robust for usage :) (do not know how slow it is, but I use it regularly and serves well)

(just replace model.tobit by generic model and data3 by name of the provided dataset and probably some functionality which would care about if there is a constatnt in the model - honestly, I am not sure if censored models without constant make any sense )

# Get coefficients
coefficients <- summary(model.tobit)$estimate[1:(length(summary(model.tobit)$estimate[,1])-1),1]
# Get names of the coefficients
names <- names(coefficients)
# Find indices of the names
indices <- unlist(lapply(names[-1], function(x) which(colnames(data3)==x)))
# Check if they exist in the dataset (if they were found)
if (length(indices)==(length(coefficients)-1)) {
# Add constant term and construct dataset
data <- cbind(rep(1, length(data3[,1])), do.call(cbind, lapply(indices, function(x) data3[,x])))
# Make the predicted values
fitted <- (data %*% coefficients)[,1]
} else {
warning("The dataset does not contain some of the variables.")
}

Best,
TK.

Followup

Message
Date: 2012-07-06 03:52
Sender: Arne Henningsen

Dear Tomas Thanks for submitting the Feature Request and for submitting code for this. I will try to implement it ASAP but I am probably not able to do this before my summer holidays. /Arne
Date: 2012-07-05 15:17
Sender: Tomas Krehlik

ok, I added truncated and censored means and put it as a function: pred <- function(model, data, type = "fitted") { coefficients <- summary(model)$estimate[1:(length(summary(model)$estimate[,1])-1),1] names <- names(coefficients) indices <- unlist(lapply(names[-1], function(x) which(colnames(data)==x))) if (length(indices)==(length(coefficients)-1)) { data_pred <- cbind(rep(1, length(data[,1])), do.call(cbind, lapply(indices, function(x) data[,x]))) fitted <- (data_pred %*% coefficients)[,1] } else { warning("The dataset does not contain some of the variables.") } sigma <- summary(model.tobit)$estimate[length(summary(model.tobit)$estimate[,1]),1] if (type == "fitted") { fitted } else if ((type == "censored") || (type == "truncated")) { if (model$left == -Inf) { print("Right censored/truncated mean:") right <- model$right truncated <- fitted - sigma*(lambda((fitted - right)/sigma)) prob_trunc <- 1 - pnorm((fitted - right)/sigma) censored <- right*dnorm((fitted - right)/sigma) + prob_trunc*truncated } else if (model$right == Inf) { print("Left censored/truncated mean:") left <- model$left truncated <- fitted + sigma*(lambda((fitted - left)/sigma)) prob_trunc <- pnorm((fitted - left)/sigma) censored <- left*dnorm((fitted - left)/sigma) + prob_trunc*truncated } else { print("Interval censored/truncated mean:") left <- model$left right <- model$right truncated <- fitted - sigma*((dnorm((left-fitted)/sigma)-dnorm((right-fitted)/sigma))/(pnorm((left-fitted)/sigma)-pnorm((right-fitted)/sigma))) prob_trunc <- pnorm((fitted - right)/sigma) - pnorm((fitted - left)/sigma) censored <- left*dnorm((fitted - left)/sigma) + prob_trunc*truncated + right*dnorm((fitted - right)/sigma) } if (type == "censored") { censored } else { truncated } } } lambda <- function(x) { lambda <- dnorm(x)/pnorm(x) lambda }

Attached Files:

Changes:

Field Old Value Date By
Operating SystemNone2012-07-06 03:49arne
ProductNone2012-07-06 03:49arne
assigned_tonone2012-07-06 03:49arne
Thanks to:
Vienna University of Economics and Business University of Wisconsin - Madison Powered By FusionForge