मजबूत रैखिक प्रतिगमन के लिए कौन सा फ़ंक्शन/पैकेज glmulti (यानी, glm की तरह व्यवहार करता है) के साथ काम करता है?

Question

पृष्ठभूमि: glmulti

glmulti साथ मल्टी-मॉडल अनुमान सामान्य रैखिक मॉडल के लिए स्वचालित मॉडल चयन के लिए एक आर समारोह/पैकेज है कि सभी संभव सामान्य रैखिक मॉडल एक आश्रित चर और भविष्यवक्ताओं का एक सेट दिया निर्माण करती है, उन्हें फिट बैठता है क्लासिक glm फ़ंक्शन के माध्यम से और फिर बहु-मॉडल अनुमान के लिए अनुमति देता है (उदाहरण के लिए, एआईसीसी, बीआईसी से प्राप्त मॉडल वजन का उपयोग करना)। glmulti किसी भी अन्य समारोह के साथ सिद्धांत में भी काम करता है जो गुणांक लौटाता है, मॉडल की लॉग-संभावना और मुफ्त पैरामीटर (और शायद अन्य जानकारी?) की संख्या उसी प्रारूप में glm करता है।

मेरा लक्ष्य: मजबूत त्रुटियों

साथ मल्टी-मॉडल अनुमान मैं प्रभाव बाहर बाहरी कारकों के कारण से बचने के लिए एक मात्रात्मक आश्रित चर के त्रुटियों की मजबूत मॉडलिंग के साथ glmulti उपयोग करना चाहते हैं।

उदाहरण के लिए, मैं मान सकता हूं कि रैखिक मॉडल में त्रुटियों को सामान्य वितरण के बजाय t distribution के रूप में वितरित किया जाता है। अपने कुर्टोसिस पैरामीटर के साथ टी वितरण में भारी पूंछ हो सकती है और इस प्रकार आउटलेटर्स (सामान्य वितरण की तुलना में) के लिए अधिक मजबूत होती है।

हालांकि, मैं टी वितरण दृष्टिकोण का उपयोग करने के लिए प्रतिबद्ध नहीं हूं। मैं किसी भी दृष्टिकोण से खुश हूं जो लॉग-संभावना को वापस देता है और इस प्रकार glmulti में मल्टीमोडेल दृष्टिकोण के साथ काम करता है। लेकिन इसका मतलब है कि, कि दुर्भाग्य से मैं आर में अच्छी तरह से ज्ञात मजबूत रैखिक मॉडल का उपयोग नहीं कर सकते हैं (जैसे, lmRobrobust या robustbase से lmrob से), क्योंकि वे लॉग-संभावना ढांचे के तहत काम करते हैं और इस प्रकार नहीं है के साथ काम नहीं कर सकते हैं glmulti।

समस्या: मैं एक मजबूत प्रतिगमन समारोह है कि साथ काम करता है नहीं मिल सकता है glmulti

आरआई के लिए ही मजबूत रेखीय प्रतीपगमन समारोह में पाया गया कि संचालित तहत लॉग-संभावना फ्रेमवर्क (heavyLm है से heavy पैकेज); यह टी वितरण के साथ त्रुटियों का मॉडल करता है। दुर्भाग्य से, heavyLmglmulti (कम से कम बॉक्स से बाहर नहीं) के साथ काम नहीं करता है, क्योंकि यह loglik (और संभवतः अन्य बातों के) के लिए कोई S3 विधि है।

इसे समझने के लिए:

summary(glm(stack.loss ~ ., data = stackloss)) 

glmulti हमारे साथ मल्टी-मॉडल निष्कर्ष:

library(glmulti) 
library(heavy) 

डाटासेट stackloss

head(stackloss) 

नियमित गाऊसी रेखीय मॉडल का उपयोग करना ing GLM के डिफ़ॉल्ट गाऊसी लिंक समारोह

stackloss.glmulti <- glmulti(stack.loss ~ ., data = stackloss, level=1, crit=bic) 
print(stackloss.glmulti) 
plot(stackloss.glmulti) 

टी वितरित त्रुटि (डिफ़ॉल्ट df = 4 है) glmulti साथ

summary(heavyLm(stack.loss ~ ., data = stackloss)) 

मल्टी मॉडल अनुमान heavyLm के रूप में बुला के साथ

लीनियर मॉडल फिटिंग फ़ंक्शन

stackloss.heavyLm.glmulti <- glmulti(stack.loss ~ ., 
data = stackloss, level=1, crit=bic, fitfunction=heavyLm) 

देता है निम्न त्रुटि:

Initialization... 
    Error in UseMethod("logLik") : 
    no applicable method for 'logLik' applied to an object of class "heavyLm". 

अगर मैं निम्नलिखित समारोह को परिभाषित,

logLik.heavyLm <- function(x){x$logLik} 

glmulti लॉग-संभावना प्राप्त कर सकते हैं, लेकिन फिर अगले त्रुटि होती है:

Initialization... 
    Error in .jcall(molly, "V", "supplyErrorDF", 
    as.integer(attr(logLik(fitfunc(as.formula(paste(y, : 
    method supplyErrorDF with signature ([I)V not found 

प्रश्न: ग्लोल्टी के साथ मजबूत रैखिक प्रतिगमन के लिए कौन सा फ़ंक्शन/पैकेज काम करता है (यानी, ग्लैम की तरह व्यवहार करता है)?

संभवतः glmulti के साथ काम करने heavyLm प्राप्त करने के लिए आगे कार्यों को परिभाषित करने के लिए एक तरीका है, लेकिन इस यात्रा पर आरंभ करने से पहले मैं पूछने के लिए कि क्या किसी को भी

• एक मजबूत रेखीय प्रतीपगमन समारोह के जानता चाहता था कि (ए) लॉग-संभावना ढांचे के तहत संचालित होता है और (बी) glm जैसा व्यवहार करता है (और इस प्रकार glmulti आउट-ऑफ-द-बॉक्स के साथ काम करेगा)।
• भारी एलएम पहले से ही glmulti के साथ काम कर रहा है।

किसी भी मदद की बहुत सराहना की जाती है!

Answer 1

यहां heavyLm का उपयोग कर उत्तर है। भले ही यह अपेक्षाकृत पुराना प्रश्न है, वही समस्या जो आपने अभी भी उल्लिखित है, तब भी भारी एलएम (यानी, त्रुटि संदेश Error in .jcall(molly, "V", "supplyErrorDF"…) का उपयोग करते समय होता है।

समस्या यह है कि glmulti को मॉडल की स्वतंत्रता की डिग्री की आवश्यकता होती है, जिसे आपको एक विशेषता के रूप में पारित करने की आवश्यकता होती है जो आपको logLik.heavyLm द्वारा दिए गए मान की विशेषता के रूप में प्रदान करने की आवश्यकता होती है; विवरण के लिए फंक्शन logLik के लिए प्रलेखन देखें।इसके अलावा, यह पता चला है कि मॉडल को फ़िट करने के लिए उपयोग किए गए डेटा बिंदुओं की संख्या को वापस करने के लिए आपको एक फ़ंक्शन प्रदान करने की भी आवश्यकता है, क्योंकि सूचना मानदंड (एआईसी, बीआईसी, ...) इस मूल्य पर भी निर्भर करता है। यह नीचे दिए गए कोड में nobs.heavyLm फ़ंक्शन द्वारा किया जाता है।

Initialization... 
TASK: Exhaustive screening of candidate set. 
Fitting... 
Completed. 

> print(stackloss.glmulti) 
glmulti.analysis 
Method: h/Fitting: glm/IC used: bic 
Level: 1/Marginality: FALSE 
From 8 models: 
Best IC: 117.892471265874 
Best model: 
[1] "stack.loss ~ 1 + Air.Flow + Water.Temp" 
Evidence weight: 0.709174196998897 
Worst IC: 162.083142797858 
2 models within 2 IC units. 
1 models to reach 95% of evidence weight. 

2 बीआईसी इकाइयों सीमा के भीतर इसलिए 2 मॉडलों के उत्पादन:

nobs.heavyLm <- function(mdl) mdl$dims[1] # the sample size (number of data points) 

logLik.heavyLm <- function(mdl) { 
    res <- mdl$logLik 
    attr(res, "nobs") <- nobs.heavyLm(mdl) # this is not really needed for 'glmulti', but is included to adhere to the format of 'logLik' 
    attr(res, "df") <- length(mdl$coefficients) + 1 + 1 # I am also considering the scale parameter that is estimated; see mdl$family 
    class(res) <- "logLik" 
    res 
} 

, जो जब कोड कि आपके द्वारा दी एक साथ रखा है, तो निम्न परिणाम पैदा करता है:

यहाँ कोड है ।

हालांकि एक महत्वपूर्ण टिप्पणी: मुझे यकीन नहीं है कि स्वतंत्रता की डिग्री के लिए ऊपर की अभिव्यक्ति सख्ती से सही है। एक मानक रैखिक मॉडल के लिए, स्वतंत्रता की डिग्री p + 1 के बराबर होगी, जहां पी मॉडल में पैरामीटर की संख्या है, और अतिरिक्त पैरामीटर (+ 1) "त्रुटि" भिन्नता है (जिसका उपयोग संभावना की गणना करने के लिए किया जाता है) । फंक्शन logLik.heavyLm ऊपर, यह मुझे स्पष्ट नहीं है कि किसी को "स्केल पैरामीटर" की गणना भी करनी चाहिए, जिसे heavyLm द्वारा आजादी की एक अतिरिक्त डिग्री के रूप में अनुमानित किया जाना चाहिए, और इसलिए p + 1 + 1, यदि मामला भी एक कार्य है तो यह मामला होगा इस पैरामीटर का। दुर्भाग्यवश, मैं इसकी पुष्टि नहीं कर सकता, क्योंकि मेरे पास इस संदर्भ तक पहुंच नहीं है कि heavyLm उद्धरण (डेम्पस्टर एट अल।, 1 9 80 का पेपर)। इस वजह से, मैं स्केल पैरामीटर की गिनती कर रहा हूं, जिससे मॉडल जटिलता के एक (थोड़ा अधिक) रूढ़िवादी अनुमान प्रदान करता है, जो "जटिल" मॉडल को दंडित करता है। छोटे अंतर नमूना मामले को छोड़कर, यह अंतर नगण्य होना चाहिए।