Ülesanded

Vajalikud paketid:

library(fmsb)
library(mlogit)

Ülesanne 1
Andmetabelist nimega neeme on Marko Neeme magistritöö (Neeme, 2012) andmed, mille uuriti 50-70-aastaste meeste suhtumist eesnäärmevähi skriiningtesti. Tunnus valmidus, näitab testis osalemise valmidust (tasemed pigem jah, ja pigem ei). Lisaks on tabelis ära toodud Suure Viisiku isiksuseomadused: neurootilisus (tunnus N), ekstravertsus (E), avatus kogemusele (O), sotsiaalsus (A) ja meelekindlus (C). Koostage binaarse logistilise regressiooni mudel, mis ennustab skriiningtestis osalemise valmidust Suure Viisiku isiksuseomaduste kaudu.

A. Esialgu tuleks saada ülevaade andmetest. Kasutage funktsiooni summary, argumendiks saab sisestada terve andmestiku või lihtsalt ühe veeru. Milline on oslejate keskmine vanus? Missugune on valmiduse tulemuste jaotus?

summary(neeme)
##      vanus             N               E               O        
##  Min.   :49.40   Min.   :1.000   Min.   :1.000   Min.   :1.000  
##  1st Qu.:56.00   1st Qu.:2.500   1st Qu.:2.750   1st Qu.:2.833  
##  Median :60.00   Median :3.000   Median :2.917   Median :3.000  
##  Mean   :61.17   Mean   :2.869   Mean   :2.909   Mean   :3.077  
##  3rd Qu.:67.00   3rd Qu.:3.250   3rd Qu.:3.167   3rd Qu.:3.333  
##  Max.   :74.48   Max.   :4.167   Max.   :4.417   Max.   :4.917  
##  NA's   :8       NA's   :3       NA's   :5       NA's   :7      
##        A               C              valmidus  
##  Min.   :1.000   Min.   :1.000   pigem ei :288  
##  1st Qu.:2.917   1st Qu.:3.083   pigem jah: 74  
##  Median :3.167   Median :3.333   NA's     :  9  
##  Mean   :3.252   Mean   :3.468                  
##  3rd Qu.:3.667   3rd Qu.:3.917                  
##  Max.   :4.750   Max.   :5.000                  
##  NA's   :5       NA's   :7
sum(is.na(neeme) == T)# kui palju puuduvaid väärtusi on kokku
## [1] 44
table(neeme$valmidus, useNA = "ifany")
## 
##  pigem ei pigem jah      <NA> 
##       288        74         9

B. Tehke ühe tunnuse tulemustest histogramm (vt. praktikum 2).

hist(neeme$N, breaks = 20, col = "gray", main = "Histogramm", ylab = "Sagedus", xlab = "Neurootlisus")

#Näide ggplot2 kasutamisest
library(ggplot2)
stacked <- stack(neeme[,2:7])
ggplot(data = stacked, aes(values))+
        geom_histogram()+
        facet_wrap(~ind)

C. Millised isiksusomadused omavad olulist seost valmidusega? Kas need omadused suurendavad või vähendavad testis osalemise valmidust?

neeme <- na.omit(neeme)# teeme uue andmetabeli, kus ei ole NA väärtusi. Kahe mudeli võrdlusega võib tekkida muidu probleeme.
neeme.m1 <- glm(valmidus ~ scale(N) + scale(E) + scale(O) + scale(A) + scale(C), data = neeme, family = binomial())
summary(neeme.m1)
## 
## Call:
## glm(formula = valmidus ~ scale(N) + scale(E) + scale(O) + scale(A) + 
##     scale(C), family = binomial(), data = neeme)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -1.6299  -0.5904  -0.4243  -0.2235   2.8783  
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.8426     0.1827 -10.088  < 2e-16 ***
## scale(N)      1.3929     0.2226   6.257 3.92e-10 ***
## scale(E)     -0.5987     0.1677  -3.570 0.000357 ***
## scale(O)      0.0680     0.1976   0.344 0.730775    
## scale(A)     -0.0340     0.1852  -0.184 0.854375    
## scale(C)      1.1171     0.2388   4.677 2.91e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 338.79  on 342  degrees of freedom
## Residual deviance: 260.14  on 337  degrees of freedom
## AIC: 272.14
## 
## Number of Fisher Scoring iterations: 5

D. Arvutage välja riskisuhted ja nende 95%-usaldusvahemikud? Milline omadus mõjutab osalemisvalmidust kõige tugevamini?

exp(coef(neeme.m1))
## (Intercept)    scale(N)    scale(E)    scale(O)    scale(A)    scale(C) 
##   0.1584011   4.0265636   0.5495177   1.0703628   0.9665714   3.0558138
exp(confint(neeme.m1))
##                 2.5 %    97.5 %
## (Intercept) 0.1085694 0.2227606
## scale(N)    2.6662030 6.3979422
## scale(E)    0.3937252 0.7632718
## scale(O)    0.7370311 1.6008183
## scale(A)    0.6755848 1.4087586
## scale(C)    1.9432587 4.9756235

E. Arvutage mudeli sobitusastet näitav hii-ruut-statistik ja selle p-väärtus.

neeme.nullmudel <- glm(valmidus~ 1, data = neeme, family = binomial())
anova(neeme.nullmudel, neeme.m1, test = "Chisq")
## Analysis of Deviance Table
## 
## Model 1: valmidus ~ 1
## Model 2: valmidus ~ scale(N) + scale(E) + scale(O) + scale(A) + scale(C)
##   Resid. Df Resid. Dev Df Deviance  Pr(>Chi)    
## 1       342     338.79                          
## 2       337     260.14  5   78.656 1.603e-15 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

F. Leidke mudeli pseudo-determinatsioonikordja.

NagelkerkeR2(neeme.m1)
## $N
## [1] 343
## 
## $R2
## [1] 0.3265317

Ülesanne 2
Tabelis ESS on lisaks tunnustele partei ja vanus veel tunnused sugu ja aastaid_koolis (vastaja kooliskäidud aastate arv). Koostage uus multinomiaalne logistiline mudel, milles on lisaks vastaja vanusele täiendavateks prediktoriteks ka sugu ja kooliskäidud aastate arv. Leidke ka riskisuhted ja nende usaldusvahemikud. Kas need muutujad seostuvad erakondlike eelistustega? Kui jah, siis millise partei puhul ja millises suunas?

ess.m2 <- mlogit(partei ~ 1| vanus + sugu + aastaid_koolis, data = ESS, shape = "wide")
summary(ess.m2)
## 
## Call:
## mlogit(formula = partei ~ 1 | vanus + sugu + aastaid_koolis, 
##     data = ESS, shape = "wide", method = "nr", print.level = 0)
## 
## Frequencies of alternatives:
##     IRL    Kesk  Reform 
## 0.19940 0.35801 0.44260 
## 
## nr method
## 4 iterations, 0h:0m:0s 
## g'(-H)^-1g = 6.7E-08 
## gradient close to zero 
## 
## Coefficients :
##                         Estimate Std. Error t-value  Pr(>|t|)    
## Kesk:(intercept)       1.7419842  0.6926221  2.5151   0.01190 *  
## Reform:(intercept)     1.5438115  0.6368543  2.4241   0.01535 *  
## Kesk:vanus             0.0114439  0.0069208  1.6536   0.09822 .  
## Reform:vanus          -0.0097754  0.0065373 -1.4953   0.13483    
## Kesk:suguN             0.3826106  0.2277829  1.6797   0.09301 .  
## Reform:suguN           0.3926685  0.2140748  1.8343   0.06662 .  
## Kesk:aastaid_koolis   -0.1538409  0.0365401 -4.2102 2.552e-05 ***
## Reform:aastaid_koolis -0.0367646  0.0332952 -1.1042   0.26951    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Log-Likelihood: -664.52
## McFadden R^2:  0.04402 
## Likelihood ratio test : chisq = 61.197 (p.value = 2.5699e-11)
# risiksuhted
data.frame(exp(coef(ess.m2)))
##                       exp.coef.ess.m2..
## Kesk:(intercept)              5.7086594
## Reform:(intercept)            4.6824035
## Kesk:vanus                    1.0115096
## Reform:vanus                  0.9902722
## Kesk:suguN                    1.4661070
## Reform:suguN                  1.4809274
## Kesk:aastaid_koolis           0.8574084
## Reform:aastaid_koolis         0.9639030
# usaldusvahemikud
exp(confint(ess.m2))
##                           2.5 %    97.5 %
## Kesk:(intercept)      1.4688355 22.186822
## Reform:(intercept)    1.3439324 16.313992
## Kesk:vanus            0.9978817  1.025324
## Reform:vanus          0.9776649  1.003042
## Kesk:suguN            0.9381591  2.291157
## Reform:suguN          0.9734485  2.252965
## Kesk:aastaid_koolis   0.7981504  0.921066
## Reform:aastaid_koolis 0.9030096  1.028903

Tabelid

library(psych)
library(htmlTable)
options(digits=2)
vtable <- describeBy(neeme)
vmat <- as.matrix(round(vtable,  digits = 2))
htmlTable(vmat, caption=" ",col.columns = c("none", "#F7F7F7"),  
          css.cell = "padding-left: .7em; padding-right: .7em;", 
          tfoot="&dagger; N - neurootilisus; E - ekstravertsus; O - avatus; A - sotsiaalsus; C - meelekindlus")
vars n mean sd median trimmed mad min max range skew kurtosis se
vanus 1 343 61.1 6.15 60 60.98 7.41 49.4 74.48 25.08 0.17 -1.1 0.33
N 2 343 2.84 0.58 3 2.87 0.49 1 4.17 3.17 -0.61 0.53 0.03
E 3 343 2.91 0.48 2.92 2.94 0.37 1 4.42 3.42 -1 3.6 0.03
O 4 343 3.08 0.47 3 3.07 0.37 1 4.92 3.92 -0.04 4.1 0.03
A 5 343 3.27 0.59 3.25 3.27 0.52 1 4.75 3.75 -0.38 1.6 0.03
C 6 343 3.49 0.63 3.33 3.48 0.62 1 5 4 -0.2 0.9 0.03
valmidus* 7 343 1.2 0.4 1 1.12 0 1 2 1 1.53 0.34 0.02
† N - neurootilisus; E - ekstravertsus; O - avatus; A - sotsiaalsus; C - meelekindlus 
library(pastecs)
options(digits=2)
stat.desc(neeme, desc=F)
##          vanus     N     E      O      A    C valmidus
## nbr.val    343 343.0 343.0  343.0  343.0  343       NA
## nbr.null     0   0.0   0.0    0.0    0.0    0       NA
## nbr.na       0   0.0   0.0    0.0    0.0    0       NA
## min         49   1.0   1.0    1.0    1.0    1       NA
## max         74   4.2   4.4    4.9    4.8    5       NA
## range       25   3.2   3.4    3.9    3.8    4       NA
## sum      20958 975.3 999.0 1057.4 1121.4 1197       NA