facto
library(dominanceanalysis)
data(longley)
head(longley)
## GNP.deflator GNP Unemployed Armed.Forces Population Year Employed
## 1947 83.0 234.289 235.6 159.0 107.608 1947 60.323
## 1948 88.5 259.426 232.5 145.6 108.632 1948 61.122
## 1949 88.2 258.054 368.2 161.6 109.773 1949 60.171
## 1950 89.5 284.599 335.1 165.0 110.929 1950 61.187
## 1951 96.2 328.975 209.9 309.9 112.075 1951 63.221
## 1952 98.1 346.999 193.2 359.4 113.270 1952 63.639
lm.1<-lm(Employed~.,longley)
da<-dominanceAnalysis(lm.1)
print(da)
##
## Dominance analysis
## Predictors: GNP.deflator, GNP, Unemployed, Armed.Forces, Population, Year
## Fit-indices: r2
##
## * Fit index: r2
## complete conditional general
## GNP.deflator Unmp,Ar.F,Pplt
## GNP Pplt GNP.,Pplt GNP.,Unmp,Ar.F,Pplt,Year
## Unemployed Ar.F
## Armed.Forces
## Population Unmp,Ar.F
## Year GNP.,Pplt GNP.,Unmp,Ar.F,Pplt
##
## Average contribution:
## GNP Year GNP.deflator Population Unemployed Armed.Forces
## 0.230 0.219 0.214 0.212 0.070 0.051
summary(da)
##
## * Fit index: r2
##
## Average contribution of each variable:
##
## GNP Year GNP.deflator Population Unemployed Armed.Forces
## 0.230 0.219 0.214 0.212 0.070 0.051
##
## Dominance Analysis matrix:
## model level fit
## 1 0 0
## GNP.deflator 1 0.943
## GNP 1 0.967
## Unemployed 1 0.253
## Armed.Forces 1 0.209
## Population 1 0.922
## Year 1 0.943
## Average level 1 1
## GNP.deflator+GNP 2 0.969
## GNP.deflator+Unemployed 2 0.959
## GNP.deflator+Armed.Forces 2 0.943
## GNP.deflator+Population 2 0.945
## GNP.deflator+Year 2 0.947
## GNP+Unemployed 2 0.981
## GNP+Armed.Forces 2 0.968
## GNP+Population 2 0.979
## GNP+Year 2 0.973
## Unemployed+Armed.Forces 2 0.561
## Unemployed+Population 2 0.969
## Unemployed+Year 2 0.982
## Armed.Forces+Population 2 0.936
## Armed.Forces+Year 2 0.947
## Population+Year 2 0.946
## Average level 2 2
## GNP.deflator+GNP+Unemployed 3 0.981
## GNP.deflator+GNP+Armed.Forces 3 0.969
## GNP.deflator+GNP+Population 3 0.982
## GNP.deflator+GNP+Year 3 0.974
## GNP.deflator+Unemployed+Armed.Forces 3 0.97
## GNP.deflator+Unemployed+Population 3 0.975
## GNP.deflator+Unemployed+Year 3 0.983
## GNP.deflator+Armed.Forces+Population 3 0.946
## GNP.deflator+Armed.Forces+Year 3 0.948
## GNP.deflator+Population+Year 3 0.948
## GNP+Unemployed+Armed.Forces 3 0.985
## GNP+Unemployed+Population 3 0.981
## GNP+Unemployed+Year 3 0.982
## GNP+Armed.Forces+Population 3 0.984
## GNP+Armed.Forces+Year 3 0.973
## GNP+Population+Year 3 0.979
## Unemployed+Armed.Forces+Population 3 0.97
## Unemployed+Armed.Forces+Year 3 0.993
## Unemployed+Population+Year 3 0.982
## Armed.Forces+Population+Year 3 0.947
## Average level 3 3
## GNP.deflator+GNP+Unemployed+Armed.Forces 4 0.985
## GNP.deflator+GNP+Unemployed+Population 4 0.982
## GNP.deflator+GNP+Unemployed+Year 4 0.983
## GNP.deflator+GNP+Armed.Forces+Population 4 0.986
## GNP.deflator+GNP+Armed.Forces+Year 4 0.974
## GNP.deflator+GNP+Population+Year 4 0.983
## GNP.deflator+Unemployed+Armed.Forces+Population 4 0.981
## GNP.deflator+Unemployed+Armed.Forces+Year 4 0.993
## GNP.deflator+Unemployed+Population+Year 4 0.983
## GNP.deflator+Armed.Forces+Population+Year 4 0.949
## GNP+Unemployed+Armed.Forces+Population 4 0.987
## GNP+Unemployed+Armed.Forces+Year 4 0.995
## GNP+Unemployed+Population+Year 4 0.983
## GNP+Armed.Forces+Population+Year 4 0.984
## Unemployed+Armed.Forces+Population+Year 4 0.995
## Average level 4 4
## GNP.deflator+GNP+Unemployed+Armed.Forces+Population 5 0.987
## GNP.deflator+GNP+Unemployed+Armed.Forces+Year 5 0.995
## GNP.deflator+GNP+Unemployed+Population+Year 5 0.984
## GNP.deflator+GNP+Armed.Forces+Population+Year 5 0.987
## GNP.deflator+Unemployed+Armed.Forces+Population+Year 5 0.995
## GNP+Unemployed+Armed.Forces+Population+Year 5 0.995
## Average level 5 5
## GNP.deflator+GNP+Unemployed+Armed.Forces+Population+Year 6 0.995
## GNP.deflator GNP Unemployed Armed.Forces Population Year
## 0.943 0.967 0.253 0.209 0.922 0.943
## 0.026 0.016 0 0.002 0.005
## 0.001 0.013 0 0.012 0.006
## 0.706 0.728 0.308 0.716 0.73
## 0.734 0.759 0.352 0.727 0.738
## 0.023 0.057 0.047 0.013 0.023
## 0.004 0.03 0.039 0.003 0.002
## 0.294 0.32 0.093 0.065 0.292 0.3
## 0.012 0.001 0.014 0.005
## 0.022 0.011 0.016 0.024
## 0.027 0.028 0.004 0.006
## 0.037 0.03 0.001 0.003
## 0.026 0.036 0.001 0.001
## 0 0.004 0.001 0.002
## 0.001 0.017 0.016 0.006
## 0.003 0.002 0.004 0
## 0 0.009 0 0.006
## 0.409 0.424 0.409 0.432
## 0.006 0.012 0.001 0.013
## 0.001 0 0.011 0
## 0.011 0.048 0.034 0.012
## 0.002 0.027 0.046 0
## 0.002 0.034 0.037 0.002
## 0.044 0.066 0.025 0.004 0.047 0.05
## 0.005 0.002 0.002
## 0.016 0.017 0.004
## 0 0.004 0
## 0.01 0 0.009
## 0.015 0.011 0.023
## 0.007 0.006 0.008
## 0 0.01 0
## 0.04 0.035 0.002
## 0.025 0.044 0
## 0.035 0.035 0.001
## 0 0.002 0.01
## 0.001 0.006 0.001
## 0.001 0.013 0
## 0.003 0.004 0
## 0 0.022 0.01
## 0.003 0.003 0.004
## 0.011 0.018 0.025
## 0 0.003 0.002
## 0.001 0 0.012
## 0.001 0.036 0.048
## 0.002 0.018 0.022 0.006 0.005 0.008
## 0.002 0.01
## 0.005 0.001
## 0.012 0.001
## 0.001 0.001
## 0.022 0.013
## 0.001 0.004
## 0.006 0.014
## 0.003 0.002
## 0.001 0.012
## 0.038 0.046
## 0 0.008
## 0 0
## 0.001 0.013
## 0.003 0.012
## 0 0.001
## 0.001 0.01 0.016 0.009 0.004 0.007
## 0.008
## 0
## 0.012
## 0.009
## 0.001
## 0
## 0 0.001 0.009 0.012 0 0.008
##
plot(da,which.graph='complete')
plot(da,which.graph='conditional')
plot(da,which.graph='general')
averageContribution(da)
##
## Average Contribution by predictor
## GNP.deflator GNP Unemployed Armed.Forces Population Year
## r2 0.214 0.23 0.07 0.051 0.212 0.219
da.no.year<-dominanceAnalysis(lm.1,constants='Year')
print(da.no.year)
##
## Dominance analysis
## Predictors: GNP.deflator, GNP, Unemployed, Armed.Forces, Population
## Constants: Year
## Fit-indices: r2
##
## * Fit index: r2
## complete conditional general
## GNP.deflator
## GNP Pplt GNP.,Pplt GNP.,Ar.F,Pplt
## Unemployed GNP,Pplt GNP.,GNP,Pplt GNP.,GNP,Ar.F,Pplt
## Armed.Forces Pplt GNP.,Pplt
## Population GNP.
##
## Average contribution:
## Unemployed GNP Armed.Forces Population GNP.deflator
## 0.025 0.016 0.007 0.002 0.001
summary(da.no.year)
##
## * Fit index: r2
##
## Average contribution of each variable:
##
## Unemployed GNP Armed.Forces Population GNP.deflator
## 0.025 0.016 0.007 0.002 0.001
##
## Dominance Analysis matrix:
## model level fit
## Year 0 0.943
## Year+GNP.deflator 1 0.947
## Year+GNP 1 0.973
## Year+Unemployed 1 0.982
## Year+Armed.Forces 1 0.947
## Year+Population 1 0.946
## Average level 1 1
## Year+GNP.deflator+GNP 2 0.974
## Year+GNP.deflator+Unemployed 2 0.983
## Year+GNP.deflator+Armed.Forces 2 0.948
## Year+GNP.deflator+Population 2 0.948
## Year+GNP+Unemployed 2 0.982
## Year+GNP+Armed.Forces 2 0.973
## Year+GNP+Population 2 0.979
## Year+Unemployed+Armed.Forces 2 0.993
## Year+Unemployed+Population 2 0.982
## Year+Armed.Forces+Population 2 0.947
## Average level 2 2
## Year+GNP.deflator+GNP+Unemployed 3 0.983
## Year+GNP.deflator+GNP+Armed.Forces 3 0.974
## Year+GNP.deflator+GNP+Population 3 0.983
## Year+GNP.deflator+Unemployed+Armed.Forces 3 0.993
## Year+GNP.deflator+Unemployed+Population 3 0.983
## Year+GNP.deflator+Armed.Forces+Population 3 0.949
## Year+GNP+Unemployed+Armed.Forces 3 0.995
## Year+GNP+Unemployed+Population 3 0.983
## Year+GNP+Armed.Forces+Population 3 0.984
## Year+Unemployed+Armed.Forces+Population 3 0.995
## Average level 3 3
## Year+GNP.deflator+GNP+Unemployed+Armed.Forces 4 0.995
## Year+GNP.deflator+GNP+Unemployed+Population 4 0.984
## Year+GNP.deflator+GNP+Armed.Forces+Population 4 0.987
## Year+GNP.deflator+Unemployed+Armed.Forces+Population 4 0.995
## Year+GNP+Unemployed+Armed.Forces+Population 4 0.995
## Average level 4 4
## Year+GNP.deflator+GNP+Unemployed+Armed.Forces+Population 5 0.995
## GNP.deflator GNP Unemployed Armed.Forces Population
## 0.004 0.03 0.039 0.003 0.002
## 0.026 0.036 0.001 0.001
## 0 0.009 0 0.006
## 0.001 0 0.011 0
## 0.002 0.027 0.046 0
## 0.002 0.034 0.037 0.002
## 0.001 0.022 0.032 0.003 0.002
## 0.01 0 0.009
## 0 0.01 0
## 0.025 0.044 0
## 0.035 0.035 0.001
## 0.001 0.013 0
## 0 0.022 0.01
## 0.003 0.003 0.004
## 0 0.003 0.002
## 0.001 0 0.012
## 0.001 0.036 0.048
## 0.001 0.017 0.027 0.007 0.004
## 0.012 0.001
## 0.022 0.013
## 0.001 0.004
## 0.003 0.002
## 0.001 0.012
## 0.038 0.046
## 0 0
## 0.001 0.013
## 0.003 0.012
## 0 0.001
## 0.001 0.011 0.02 0.01 0.004
## 0
## 0.012
## 0.009
## 0.001
## 0
## 0 0.001 0.009 0.012 0
##
plot(da.no.year,which.graph='complete')
plot(da.no.year,which.graph='conditional')
da.terms=c(GNP.rel='GNP.deflator+GNP',
pop.rel='Unemployed+Armed.Forces+Population+Unemployed',
year='Year')
da.grouped<-dominanceAnalysis(lm.1,terms=da.terms)
print(da.grouped)
##
## Dominance analysis
## Predictors: GNP.deflator+GNP, Unemployed+Armed.Forces+Population+Unemployed, Year
## Terms: GNP.rel = GNP.deflator+GNP ; pop.rel = Unemployed+Armed.Forces+Population+Unemployed ; year = Year
## Fit-indices: r2
##
## * Fit index: r2
## complete conditional general
## GNP.rel year
## pop.rel GNP.,year GNP.,year GNP.,year
## year
##
## Average contribution:
## pop.rel GNP.rel year
## 0.342 0.331 0.322
summary(da.grouped)
##
## * Fit index: r2
##
## Average contribution of each variable:
##
## pop.rel GNP.rel year
## 0.342 0.331 0.322
##
## Dominance Analysis matrix:
## model level fit GNP.rel pop.rel year
## 1 0 0 0.969 0.97 0.943
## GNP.rel 1 0.969 0.019 0.005
## pop.rel 1 0.97 0.018 0.025
## year 1 0.943 0.03 0.051
## Average level 1 1 0.024 0.035 0.015
## GNP.rel+pop.rel 2 0.987 0.008
## GNP.rel+year 2 0.974 0.022
## pop.rel+year 2 0.995 0.001
## Average level 2 2 0.001 0.022 0.008
## GNP.rel+pop.rel+year 3 0.995
plot(da.grouped, which.graph='conditional')
data("tropicbird")
head(tropicbird)
## ID rem alt slo rain coast land pres
## 1 1 5.76706 67 9.99194 6395.874 0.00657 2 0
## 2 2 5.77885 49 10.27684 6397.938 0.00471 2 0
## 3 3 5.79152 32 7.63426 6400.896 0.00287 2 0
## 4 4 5.80387 23 7.42965 6403.822 0.00108 2 0
## 5 5 5.82159 53 12.27207 6455.543 0.00231 2 0
## 6 6 5.26972 243 11.22250 6901.880 0.00443 2 0
library(caTools)
set.seed(101)
sample <- caTools::sample.split(tropicbird$ID, SplitRatio = .70)
train <- subset(tropicbird, sample == TRUE)
test <- subset(tropicbird, sample == FALSE)
modpres <- glm(pres~rem+land+alt+slo+rain+coast, data=train, family=binomial(link='logit'))
summary(modpres)
##
## Call:
## glm(formula = pres ~ rem + land + alt + slo + rain + coast, family = binomial(link = "logit"),
## data = train)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.5118 -0.2137 -0.0901 -0.0298 3.2988
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.255e+01 1.669e+00 -7.518 5.54e-14 ***
## rem 6.225e-01 1.998e-01 3.116 0.001833 **
## land 3.009e-01 4.048e-01 0.743 0.457210
## alt 1.920e-03 4.835e-04 3.972 7.13e-05 ***
## slo 1.588e-02 1.567e-02 1.013 0.310834
## rain 6.892e-04 1.706e-04 4.040 5.35e-05 ***
## coast 2.865e+01 8.476e+00 3.380 0.000725 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 530.66 on 1677 degrees of freedom
## Residual deviance: 352.15 on 1671 degrees of freedom
## AIC: 366.15
##
## Number of Fisher Scoring iterations: 8
anova(modpres, test="Chisq")
## Analysis of Deviance Table
##
## Model: binomial, link: logit
##
## Response: pres
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 1677 530.66
## rem 1 145.298 1676 385.36 < 2.2e-16 ***
## land 1 1.446 1675 383.92 0.2292467
## alt 1 10.647 1674 373.27 0.0011027 **
## slo 1 0.718 1673 372.55 0.3967159
## rain 1 8.057 1672 364.49 0.0045317 **
## coast 1 12.346 1671 352.15 0.0004418 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
library(pscl)
## Classes and Methods for R developed in the
## Political Science Computational Laboratory
## Department of Political Science
## Stanford University
## Simon Jackman
## hurdle and zeroinfl functions by Achim Zeileis
pR2(modpres)
## fitting null model for pseudo-r2
## llh llhNull G2 McFadden r2ML r2CU
## -176.0740171 -265.3300722 178.5121103 0.3363963 0.1009205 0.3722362
dapres<-dominanceAnalysis(modpres)
dominanceMatrix(dapres, type="complete",
fit.functions = "r2.m", ordered=TRUE)
## rem coast alt rain land slo
## rem 0.5 0.5 0.5 0.5 1.0 1.0
## coast 0.5 0.5 0.5 0.5 0.5 1.0
## alt 0.5 0.5 0.5 0.5 0.5 0.5
## rain 0.5 0.5 0.5 0.5 0.5 0.5
## land 0.0 0.5 0.5 0.5 0.5 0.5
## slo 0.0 0.0 0.5 0.5 0.5 0.5
plot(dapres, which.graph ="conditional",
fit.function = "r2.m")
