trviews <-
function(x, Y, ii = 1)
{
# Trimmed views for 90, 80, ... 0 percent
# trimming. Increase ii if 90% trimming is too harsh.
# Allows visualization of m and crudely estimation of
# c beta in models of the form y = m(x^T beta) + e.
# Workstation: activate a graphics device
# with commands "X11()" or "motif()."
# R needs command "library(lqs)."
# Advance the view with the right mouse button and
# in R, highight "stop."
x <- as.matrix(x)
out <- cov.mcd(x)
center <- out$center
cov <- out$cov
rd2 <- mahalanobis(x, center, cov)
labs <- c("90%", "80%", "70%", "60%", "50%", "40%", "30%",
"20%","10%","0%")
tem <- seq(0.1, 1, 0.1)
for(i in ii:10) {
val <- quantile(rd2, tem[i])
b <- lsfit(x[rd2 <= val, ], Y[rd2 <= val])$coef
ESP <- x %*% b[-1]
plot(ESP, Y)
title(labs[i])
identify(ESP, Y)
print(b)
}
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.