predict.lm में क्लस्टर सहप्रसरण मैट्रिक्स()

Question

मैं एक डाटासेट जिसमें डेटा कई समूहों (क्षेत्रों में शहरों) में क्लस्टर है का विश्लेषण करने कर रहा हूँ का उपयोग करना। डाटासेट लगता है:

R> df <- data.frame(x = rnorm(10), 
        y = 3*rnorm(x), 
        groups = factor(sample(c('0','1'), 10, TRUE))) 
R> head(df) 
     x  y groups 
1 -0.8959 1.54  1 
2 -0.1008 -2.73  1 
3 0.4406 0.44  0 
4 0.0683 1.62  1 
5 -0.0037 -0.20  1 
6 -0.8966 -2.34  0 

मैं चाहता हूँ मेरी एलएम() समूहों में intraclass सहसंबंध के लिए खाते का अनुमान है और उस उद्देश्य के लिए मैं एक समारोह cl() कि एक lm() लेता है और मजबूत क्लस्टर सहप्रसरण मैट्रिक्स रिटर्न (मूल here उपयोग कर रहा हूँ):

cl <- function(fm, cluster) { 
    library(sandwich) 
    M <- length(unique(cluster)) 
    N <- length(cluster)    
    K <- fm$rank     
    dfc <- (M/(M-1))*((N-1)/(N-K-1)) 
    uj <- apply(estfun(fm), 2, function(x) tapply(x, cluster, sum)); 
    vcovCL <- dfc * sandwich(fm, meat = crossprod(uj)/N) 
    return(vcovCL) 
} 

अब,

output <- lm(y ~ x, data = df) 
clcov <- cl(output, df$groups) 
coeftest(output, clcov, nrow(df) - 1) 

मुझे अनुमान मैं की जरूरत है देता है। समस्या अब है कि मैं भविष्यवाणी के लिए मॉडल का उपयोग करना चाहते हैं, और मैं भविष्यवाणी की मानक त्रुटि नई सहप्रसरण मैट्रिक्स clcov के साथ गणना की जानी चाहिए। जो है, मैं

predict(output, se.fit = TRUE) 

की जरूरत है लेकिन vcov(output) के बजाय clcov का उपयोग कर। vcov() <- की तरह कुछ सही होगा।

बेशक, मैं भविष्यवाणियों के लिए अपना खुद का कार्य लिख सकता हूं, लेकिन मैं सोच रहा हूं कि एक और व्यावहारिक विधि है जो मुझे lm (जैसे arm :: sim) के लिए विधियों का उपयोग करने की अनुमति देती है।

Answer 1

भविष्यवाणी में se.fit vcov मैट्रिक्स का उपयोग करते हुए, लेकिन qr अपघटन और अवशिष्ट विचरण का उपयोग कर की गणना नहीं की है। यह रूप में अच्छी तरह vcov() फ़ंक्शन के लिए चला जाता है: यह() एक साथ अवशिष्ट विचरण के साथ summary.lm से बगैर माप cov मैट्रिक्स लेता है, और उन लोगों का उपयोग करता है। और unscaled cov matrix है - फिर से- क्यूआर अपघटन से गणना की।

तो मुझे डर है कि उत्तर "नहीं, वहाँ कोई अन्य विकल्प अपने खुद के समारोह लिखने के लिए की तुलना में है" कर रहा हूँ। आप वास्तव में vcov मैट्रिक्स सेट नहीं कर सकते हैं, क्योंकि इसकी आवश्यकता होने पर इसे फिर से गणना की जाती है। फिर भी, अपना खुद का कार्य लिखना तुच्छ है।

predict.rob <- function(x,clcov,newdata){ 
    if(missing(newdata)){ newdata <- x$model } 
    m.mat <- model.matrix(x$terms,data=newdata) 
    m.coef <- x$coef 
    fit <- as.vector(m.mat %*% x$coef) 
    se.fit <- sqrt(diag(m.mat%*%clcov%*%t(m.mat))) 
    return(list(fit=fit,se.fit=se.fit)) 
} 

मुझे लगता है() फ़ंक्शन का उपयोग नहीं किया अनावश्यक गणना से बचने के लिए। यह वैसे भी कोड को बहुत छोटा नहीं करेगा।

---

एक तरफ ध्यान दें पर, इस तरह सवाल बेहतर stats.stackexchange.com

Answer 2

पर कहा जाता है मैं थोड़ा भविष्यवाणी समारोह के साथ और स्थायी उपरोक्त कोड संशोधित - इस तरह से आप के लिए मान दर्ज करने की उम्मीद नहीं कर रहे हैं न्यूडाटा डेटाफ्रेम

predict.rob <- function(x,clcov,newdata){ 
if(missing(newdata)){ newdata <- x$model } 
tt <- terms(x) 
Terms <- delete.response(tt) 
m.mat <- model.matrix(Terms,data=newdata) 
m.coef <- x$coef 
fit <- as.vector(m.mat %*% x$coef) 
se.fit <- sqrt(diag(m.mat%*%clcov%*%t(m.mat))) 
return(list(fit=fit,se.fit=se.fit))}