1 Session 05: Binary logistic regression modeling 2

Does the choice of a genitive construction (of vs. s) vary as a function of

  • a length-based short-before-long effect, but, this time, if we hypothesize a short-before-long effect, maybe we should not just be looking at the length of the possessor, but how the length of the possessor compares to the length of the possessum (which is a variable we may have to create from the data);
  • the degree/kind of animacy of the possessor (POR_ANIMACY: animate vs. collective vs. inanimate vs. locative vs. temporal);
  • whether the speakers are non-native speakers or native speakers of English (SPEAKER: nns vs. ns);
  • any pairwise interaction of these predictors;
  • the three-way interaction of these predictors?
rm(list=ls(all.names=TRUE))
library(car); library(effects)
source("105_04-05_R2.r"); source("105_04-05_C.score.r") # load small helper functions
summary(x <- read.delim(      # summarize x, the result of loading
   file="105_04-05_GENs.csv", # this file
   stringsAsFactors=TRUE))    # change categorical variables into factors
##       CASE      SPEAKER       MODALITY    GENITIVE    POR_LENGTH    
##  Min.   :   2   nns:2666   spoken :1685   of:2720   Min.   :  1.00  
##  1st Qu.:1006   ns : 934   written:1915   s : 880   1st Qu.:  8.00  
##  Median :2018                                       Median : 11.00  
##  Mean   :2012                                       Mean   : 14.58  
##  3rd Qu.:3017                                       3rd Qu.: 17.00  
##  Max.   :4040                                       Max.   :204.00  
##    PUM_LENGTH         POR_ANIMACY   POR_FINAL_SIB        POR_DEF    
##  Min.   :  2.00   animate   : 920   absent :2721   definite  :2349  
##  1st Qu.:  6.00   collective: 607   present: 879   indefinite:1251  
##  Median :  9.00   inanimate :1671                                   
##  Mean   : 10.35   locative  : 243                                   
##  3rd Qu.: 13.00   temporal  : 159                                   
##  Max.   :109.00

Let’s already compute the baselines for what will be the response variable, GENITIVE:

(baseline.1 <- max(        # make baseline.1 the highest
   prop.table(             # proportion in the
      table(x$GENITIVE)))) # frequency table of the response variable
## [1] 0.7555556
(baseline.2 <- sum(          # make baseline.2 the sum of the
   prop.table(               # proportions in the
      table(x$GENITIVE))^2)) # frequency table of the response variable squared
## [1] 0.6306173

1.1 Deviance

Let’s compute the deviance of the null model:

m.00 <- glm(GENITIVE ~ 1, data=x, family=binomial, na.action=na.exclude)
deviance(m.00)
## [1] 4004.273
sum(residuals(m.00)^2)
## [1] 4004.273

1.2 Exploration & preparation

Some exploration of the relevant variables:

# the predictor(s)/response on its/their own
hist(x$LEN_PORmPUM_LOG <- log2(x$POR_LENGTH)-log2(x$PUM_LENGTH), main="")

table(x$POR_ANIMACY)
## 
##    animate collective  inanimate   locative   temporal 
##        920        607       1671        243        159
table(x$SPEAKER)
## 
##  nns   ns 
## 2666  934
# the predictor(s) w/ the response
spineplot(x$GENITIVE ~ x$LEN_PORmPUM_LOG)

table(x$POR_ANIMACY, x$GENITIVE)
##             
##                of    s
##   animate     370  550
##   collective  408  199
##   inanimate  1638   33
##   locative    199   44
##   temporal    105   54
table(x$SPEAKER, x$GENITIVE)
##      
##         of    s
##   nns 2024  642
##   ns   696  238
ftable(x$POR_ANIMACY, x$SPEAKER, x$GENITIVE)
##                   of    s
##                          
## animate    nns   272  406
##            ns     98  144
## collective nns   314  138
##            ns     94   61
## inanimate  nns  1250   27
##            ns    388    6
## locative   nns   117   34
##            ns     82   10
## temporal   nns    71   37
##            ns     34   17

1.3 A multifactorial model selection process

Let’s fit our initial regression model:

summary(m.01 <- glm(      # make/summarize the gen. linear model m.01:
   GENITIVE ~ 1 +         # GENITIVE ~ an overall intercept (1)
   LEN_PORmPUM_LOG*POR_ANIMACY*SPEAKER, # & these predictors & their interaction
   data=x,                # those variables are in x
   family=binomial,       # the response is binary
   na.action=na.exclude)) # skip cases with NA/missing data
## 
## Call:
## glm(formula = GENITIVE ~ 1 + LEN_PORmPUM_LOG * POR_ANIMACY * 
##     SPEAKER, family = binomial, data = x, na.action = na.exclude)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -2.1516  -0.4186  -0.1923  -0.0599   3.0249  
## 
## Coefficients:
##                                                 Estimate Std. Error z value
## (Intercept)                                      0.65632    0.09135   7.184
## LEN_PORmPUM_LOG                                 -0.70601    0.08154  -8.659
## POR_ANIMACYcollective                           -1.77888    0.16150 -11.015
## POR_ANIMACYinanimate                            -4.32698    0.21686 -19.953
## POR_ANIMACYlocative                             -2.30616    0.27227  -8.470
## POR_ANIMACYtemporal                             -1.32138    0.22536  -5.863
## SPEAKERns                                       -0.08278    0.17204  -0.481
## LEN_PORmPUM_LOG:POR_ANIMACYcollective           -0.46680    0.15054  -3.101
## LEN_PORmPUM_LOG:POR_ANIMACYinanimate             0.32099    0.17729   1.811
## LEN_PORmPUM_LOG:POR_ANIMACYlocative             -0.29249    0.25885  -1.130
## LEN_PORmPUM_LOG:POR_ANIMACYtemporal              0.44522    0.18375   2.423
## LEN_PORmPUM_LOG:SPEAKERns                        0.14245    0.16001   0.890
## POR_ANIMACYcollective:SPEAKERns                  0.79493    0.29353   2.708
## POR_ANIMACYinanimate:SPEAKERns                  -0.28854    0.51213  -0.563
## POR_ANIMACYlocative:SPEAKERns                   -0.50011    0.49785  -1.005
## POR_ANIMACYtemporal:SPEAKERns                    0.16161    0.41275   0.392
## LEN_PORmPUM_LOG:POR_ANIMACYcollective:SPEAKERns -0.14283    0.28943  -0.493
## LEN_PORmPUM_LOG:POR_ANIMACYinanimate:SPEAKERns  -0.65372    0.43153  -1.515
## LEN_PORmPUM_LOG:POR_ANIMACYlocative:SPEAKERns   -0.19485    0.50679  -0.384
## LEN_PORmPUM_LOG:POR_ANIMACYtemporal:SPEAKERns   -0.53061    0.38091  -1.393
##                                                 Pr(>|z|)    
## (Intercept)                                     6.75e-13 ***
## LEN_PORmPUM_LOG                                  < 2e-16 ***
## POR_ANIMACYcollective                            < 2e-16 ***
## POR_ANIMACYinanimate                             < 2e-16 ***
## POR_ANIMACYlocative                              < 2e-16 ***
## POR_ANIMACYtemporal                             4.54e-09 ***
## SPEAKERns                                        0.63039    
## LEN_PORmPUM_LOG:POR_ANIMACYcollective            0.00193 ** 
## LEN_PORmPUM_LOG:POR_ANIMACYinanimate             0.07022 .  
## LEN_PORmPUM_LOG:POR_ANIMACYlocative              0.25849    
## LEN_PORmPUM_LOG:POR_ANIMACYtemporal              0.01539 *  
## LEN_PORmPUM_LOG:SPEAKERns                        0.37331    
## POR_ANIMACYcollective:SPEAKERns                  0.00677 ** 
## POR_ANIMACYinanimate:SPEAKERns                   0.57315    
## POR_ANIMACYlocative:SPEAKERns                    0.31511    
## POR_ANIMACYtemporal:SPEAKERns                    0.69538    
## LEN_PORmPUM_LOG:POR_ANIMACYcollective:SPEAKERns  0.62166    
## LEN_PORmPUM_LOG:POR_ANIMACYinanimate:SPEAKERns   0.12980    
## LEN_PORmPUM_LOG:POR_ANIMACYlocative:SPEAKERns    0.70062    
## LEN_PORmPUM_LOG:POR_ANIMACYtemporal:SPEAKERns    0.16362    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 4004.3  on 3599  degrees of freedom
## Residual deviance: 2404.3  on 3580  degrees of freedom
## AIC: 2444.3
## 
## Number of Fisher Scoring iterations: 7
pchisq(              # compute the area under the chi-squared curve
   q=4004.3-2404.3,  # for this chi-squared value: null - residual deviance
   df=3599-3580,     # at this df: null - residual df
   lower.tail=FALSE) # only using the right/upper tail/side
## [1] 0
drop1(           # drop each droppable predictor at a time
   m.01,         # from the model m.02 &
   test="Chisq") # test its significance w/ a LRT
## Single term deletions
## 
## Model:
## GENITIVE ~ 1 + LEN_PORmPUM_LOG * POR_ANIMACY * SPEAKER
##                                     Df Deviance    AIC    LRT Pr(>Chi)
## <none>                                   2404.3 2444.3                
## LEN_PORmPUM_LOG:POR_ANIMACY:SPEAKER  4   2408.0 2440.0 3.7384   0.4426

Let’s fit the 2nd model with that 3-way interaction LEN_PORmPUM_LOG:POR_ANIMACY:SPEAKER deleted:

summary(m.02 <- update( # make m.02 an update of
   m.01, .~.            # m.01, namely all of it (.~.), but then
   - LEN_PORmPUM_LOG:POR_ANIMACY:SPEAKER)) # remove this interaction
## 
## Call:
## glm(formula = GENITIVE ~ LEN_PORmPUM_LOG + POR_ANIMACY + SPEAKER + 
##     LEN_PORmPUM_LOG:POR_ANIMACY + LEN_PORmPUM_LOG:SPEAKER + POR_ANIMACY:SPEAKER, 
##     family = binomial, data = x, na.action = na.exclude)
## 
## Deviance Residuals: 
##      Min        1Q    Median        3Q       Max  
## -2.10679  -0.41718  -0.18926  -0.09363   3.07952  
## 
## Coefficients:
##                                       Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                            0.63865    0.08956   7.131 9.99e-13 ***
## LEN_PORmPUM_LOG                       -0.66568    0.07506  -8.869  < 2e-16 ***
## POR_ANIMACYcollective                 -1.75874    0.15934 -11.038  < 2e-16 ***
## POR_ANIMACYinanimate                  -4.30116    0.21500 -20.006  < 2e-16 ***
## POR_ANIMACYlocative                   -2.29342    0.26298  -8.721  < 2e-16 ***
## POR_ANIMACYtemporal                   -1.31543    0.22698  -5.795 6.82e-09 ***
## SPEAKERns                             -0.01779    0.16826  -0.106   0.9158    
## LEN_PORmPUM_LOG:POR_ANIMACYcollective -0.50145    0.12853  -3.901 9.56e-05 ***
## LEN_PORmPUM_LOG:POR_ANIMACYinanimate   0.20466    0.15951   1.283   0.1995    
## LEN_PORmPUM_LOG:POR_ANIMACYlocative   -0.34103    0.22263  -1.532   0.1256    
## LEN_PORmPUM_LOG:POR_ANIMACYtemporal    0.31318    0.16011   1.956   0.0505 .  
## LEN_PORmPUM_LOG:SPEAKERns             -0.02120    0.11639  -0.182   0.8555    
## POR_ANIMACYcollective:SPEAKERns        0.72665    0.29266   2.483   0.0130 *  
## POR_ANIMACYinanimate:SPEAKERns        -0.30015    0.48555  -0.618   0.5365    
## POR_ANIMACYlocative:SPEAKERns         -0.55095    0.45355  -1.215   0.2245    
## POR_ANIMACYtemporal:SPEAKERns          0.08259    0.40494   0.204   0.8384    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 4004.3  on 3599  degrees of freedom
## Residual deviance: 2408.0  on 3584  degrees of freedom
## AIC: 2440
## 
## Number of Fisher Scoring iterations: 7
anova(           # compare
   m.01,         # m.01
   m.02,         # to m.02
   test="Chisq") # using a LRT
## Analysis of Deviance Table
## 
## Model 1: GENITIVE ~ 1 + LEN_PORmPUM_LOG * POR_ANIMACY * SPEAKER
## Model 2: GENITIVE ~ LEN_PORmPUM_LOG + POR_ANIMACY + SPEAKER + LEN_PORmPUM_LOG:POR_ANIMACY + 
##     LEN_PORmPUM_LOG:SPEAKER + POR_ANIMACY:SPEAKER
##   Resid. Df Resid. Dev Df Deviance Pr(>Chi)
## 1      3580     2404.3                     
## 2      3584     2408.0 -4  -3.7384   0.4426
drop1(           # drop each droppable predictor at a time
   m.02,         # from the model m.02 &
   test="Chisq") # test its significance w/ a LRT
## Single term deletions
## 
## Model:
## GENITIVE ~ LEN_PORmPUM_LOG + POR_ANIMACY + SPEAKER + LEN_PORmPUM_LOG:POR_ANIMACY + 
##     LEN_PORmPUM_LOG:SPEAKER + POR_ANIMACY:SPEAKER
##                             Df Deviance    AIC    LRT  Pr(>Chi)    
## <none>                           2408.0 2440.0                     
## LEN_PORmPUM_LOG:POR_ANIMACY  4   2439.8 2463.8 31.766 2.136e-06 ***
## LEN_PORmPUM_LOG:SPEAKER      1   2408.0 2438.0  0.033   0.85530    
## POR_ANIMACY:SPEAKER          4   2418.9 2442.9 10.880   0.02794 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Let’s fit the 3rd model with that interaction LEN_PORmPUM_LOG:SPEAKER deleted:

summary(m.03 <- update( # make m.03 an update of
   m.02, .~.            # m.02, namely all of it (.~.), but then
   - LEN_PORmPUM_LOG:SPEAKER)) # remove this interaction
## 
## Call:
## glm(formula = GENITIVE ~ LEN_PORmPUM_LOG + POR_ANIMACY + SPEAKER + 
##     LEN_PORmPUM_LOG:POR_ANIMACY + POR_ANIMACY:SPEAKER, family = binomial, 
##     data = x, na.action = na.exclude)
## 
## Deviance Residuals: 
##      Min        1Q    Median        3Q       Max  
## -2.11235  -0.41886  -0.18961  -0.09468   3.08190  
## 
## Coefficients:
##                                       Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                            0.64082    0.08887   7.211 5.57e-13 ***
## LEN_PORmPUM_LOG                       -0.67066    0.06999  -9.582  < 2e-16 ***
## POR_ANIMACYcollective                 -1.76343    0.15747 -11.198  < 2e-16 ***
## POR_ANIMACYinanimate                  -4.30316    0.21475 -20.038  < 2e-16 ***
## POR_ANIMACYlocative                   -2.29898    0.26155  -8.790  < 2e-16 ***
## POR_ANIMACYtemporal                   -1.31845    0.22659  -5.819 5.93e-09 ***
## SPEAKERns                             -0.02626    0.16145  -0.163  0.87078    
## LEN_PORmPUM_LOG:POR_ANIMACYcollective -0.50226    0.12847  -3.910 9.25e-05 ***
## LEN_PORmPUM_LOG:POR_ANIMACYinanimate   0.20638    0.15920   1.296  0.19486    
## LEN_PORmPUM_LOG:POR_ANIMACYlocative   -0.34169    0.22268  -1.534  0.12492    
## LEN_PORmPUM_LOG:POR_ANIMACYtemporal    0.31237    0.16007   1.951  0.05100 .  
## POR_ANIMACYcollective:SPEAKERns        0.73847    0.28489   2.592  0.00954 ** 
## POR_ANIMACYinanimate:SPEAKERns        -0.29379    0.48413  -0.607  0.54396    
## POR_ANIMACYlocative:SPEAKERns         -0.53295    0.44234  -1.205  0.22827    
## POR_ANIMACYtemporal:SPEAKERns          0.08969    0.40281   0.223  0.82380    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 4004.3  on 3599  degrees of freedom
## Residual deviance: 2408.0  on 3585  degrees of freedom
## AIC: 2438
## 
## Number of Fisher Scoring iterations: 7
anova(           # compare
   m.02,         # m.02
   m.03,         # to m.03
   test="Chisq") # using a LRT
## Analysis of Deviance Table
## 
## Model 1: GENITIVE ~ LEN_PORmPUM_LOG + POR_ANIMACY + SPEAKER + LEN_PORmPUM_LOG:POR_ANIMACY + 
##     LEN_PORmPUM_LOG:SPEAKER + POR_ANIMACY:SPEAKER
## Model 2: GENITIVE ~ LEN_PORmPUM_LOG + POR_ANIMACY + SPEAKER + LEN_PORmPUM_LOG:POR_ANIMACY + 
##     POR_ANIMACY:SPEAKER
##   Resid. Df Resid. Dev Df  Deviance Pr(>Chi)
## 1      3584       2408                      
## 2      3585       2408 -1 -0.033257   0.8553
drop1(           # drop each droppable predictor at a time
   m.03,         # from the model m.03 &
   test="Chisq") # test its significance w/ a LRT
## Single term deletions
## 
## Model:
## GENITIVE ~ LEN_PORmPUM_LOG + POR_ANIMACY + SPEAKER + LEN_PORmPUM_LOG:POR_ANIMACY + 
##     POR_ANIMACY:SPEAKER
##                             Df Deviance    AIC    LRT  Pr(>Chi)    
## <none>                           2408.0 2438.0                     
## LEN_PORmPUM_LOG:POR_ANIMACY  4   2440.0 2462.0 31.930 1.977e-06 ***
## POR_ANIMACY:SPEAKER          4   2419.1 2441.1 11.066   0.02583 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
m.final <- m.03; rm(m.03); invisible(gc())

Let’s make ‘predictions’ and evaluate them:

x$PREDS.NUM.2 <- predict( # make x$PREDS.NUM.2 the result of predicting
   m.final,               # from m.final
   type="response")       # predicted probabilities of "s"
# if you don't say type="response", predict gives you log odds
x$PREDS.CAT <- factor(        # make x$PREDS.CAT the result of checking
   ifelse(x$PREDS.NUM.2>=0.5, # if the predicted prob. of "s" is >=0.5
      levels(x$GENITIVE)[2],  # if yes, predict "s"
      levels(x$GENITIVE)[1])) # otherwise, predict "of"
(qwe <- table(           # confusion matrix: cross-tabulate
   "OBS"=x$GENITIVE,     # observed genitives in the rows
   "PREDS"=x$PREDS.CAT)) # predicted orders in the columns
##     PREDS
## OBS    of    s
##   of 2454  266
##   s   310  570
c(rev(tapply(x$GENITIVE==x$PREDS.CAT, x$GENITIVE, mean)), # accuracies for ...
  "overall"=mean(x$GENITIVE==x$PREDS.CAT)) # ... each level & overall
##         s        of   overall 
## 0.6477273 0.9022059 0.8400000
c( # evaluate the confusion matrix
   "Prec. for s" =qwe[ "s", "s"] / sum(qwe[    , "s"]),
   "Rec. for s"  =qwe[ "s", "s"] / sum(qwe[ "s",    ]),
   "Prec. for of"=qwe["of","of"] / sum(qwe[    ,"of"]),
   "Rec. for of" =qwe["of","of"] / sum(qwe["of",    ]))
##  Prec. for s   Rec. for s Prec. for of  Rec. for of 
##    0.6818182    0.6477273    0.8878437    0.9022059
c("R2"=R2(m.final), "C"=C.score(m.final)) # R-squared & C
##        R2         C 
## 0.5335965 0.8986036

Let’s visualize the nature of the effect of LEN_PORmPUM_LOG:POR_ANIMACY based on the predictions:

(ph.la <- data.frame(                     # make ph.la a data frame of
   len.ani <- effect(                     # the effect
      "LEN_PORmPUM_LOG:POR_ANIMACY",      # of LEN_PORmPUM_LOG:POR_ANIMACY
      xlevels=list(LEN_PORmPUM_LOG=-4:4), # w/ these values for LEN_PORmPUM_LOG
      m.final)))                          # in m.final
##    LEN_PORmPUM_LOG POR_ANIMACY         fit          se        lower       upper
## 1               -4     animate 0.964995167 0.010732807 0.9366654268 0.980911180
## 2               -3     animate 0.933762692 0.015464168 0.8962219877 0.958354439
## 3               -2     animate 0.878181385 0.019688827 0.8340490730 0.911818185
## 4               -1     animate 0.786618460 0.020577239 0.7435273242 0.824181653
## 5                0     animate 0.653396580 0.017519047 0.6183110604 0.686890163
## 6                1     animate 0.490837160 0.020513568 0.4507767242 0.531015622
## 7                2     animate 0.330192632 0.029080480 0.2758733167 0.389455090
## 8                3     animate 0.201334835 0.031196146 0.1470135317 0.269388791
## 9                4     animate 0.114190794 0.026335613 0.0718292027 0.176777133
## 10              -4  collective 0.977108479 0.009125137 0.9504633183 0.989578685
## 11              -3  collective 0.929623169 0.019952030 0.8790212514 0.960022358
## 12              -2  collective 0.803450188 0.032664294 0.7315682184 0.859772887
## 13              -1  collective 0.558501674 0.030970559 0.4972271761 0.618044933
## 14               0  collective 0.281338372 0.022128739 0.2400611313 0.326662174
## 15               1  collective 0.108056512 0.017076640 0.0788544432 0.146354725
## 16               2  collective 0.036135870 0.009475954 0.0215225418 0.060062271
## 17               3  collective 0.011468910 0.004239852 0.0055436094 0.023577327
## 18               4  collective 0.003577539 0.001705505 0.0014037321 0.009087051
## 19              -4   inanimate 0.131443260 0.069510582 0.0438729443 0.332938978
## 20              -3   inanimate 0.086864232 0.037608215 0.0361996368 0.194153942
## 21              -2   inanimate 0.056422056 0.018427345 0.0294486514 0.105418163
## 22              -1   inanimate 0.036225313 0.008239215 0.0231206688 0.056329289
## 23               0   inanimate 0.023081315 0.004063324 0.0163254197 0.032540495
## 24               1   inanimate 0.014634072 0.003229110 0.0094845324 0.022515951
## 25               2   inanimate 0.009249061 0.003021352 0.0048681010 0.017503241
## 26               3   inanimate 0.005833881 0.002646221 0.0023938177 0.014147422
## 27               4   inanimate 0.003675066 0.002161627 0.0011583720 0.011596083
## 28              -4    locative 0.904318818 0.065050659 0.6840949593 0.976331696
## 29              -3    locative 0.774484239 0.095877338 0.5393933482 0.909678349
## 30              -2    locative 0.555138390 0.087955205 0.3830574676 0.714940059
## 31              -1    locative 0.311975694 0.043846419 0.2330327080 0.403589135
## 32               0    locative 0.141455667 0.026043363 0.0976551339 0.200535737
## 33               1    locative 0.056486712 0.019914540 0.0279778871 0.110735922
## 34               2    locative 0.021290829 0.011835875 0.0070950857 0.062112445
## 35               3    locative 0.007842601 0.006002145 0.0017399124 0.034608204
## 36               4    locative 0.002864010 0.002793417 0.0004221157 0.019161222
## 37              -4    temporal 0.683948493 0.127825282 0.4044341600 0.873357610
## 38              -3    temporal 0.601972295 0.109156405 0.3824363698 0.786943563
## 39              -2    temporal 0.513847009 0.081662166 0.3577119523 0.667325461
## 40              -1    temporal 0.424852275 0.053325703 0.3250602122 0.531170088
## 41               0    temporal 0.340476161 0.038786244 0.2689980630 0.420033935
## 42               1    temporal 0.265132480 0.045062208 0.1865229243 0.362123903
## 43               2    temporal 0.201370889 0.055382698 0.1137786495 0.331195140
## 44               3    temporal 0.149817159 0.060437310 0.0650079736 0.308734224
## 45               4    temporal 0.109649779 0.059635675 0.0358613321 0.289652708
plot(len.ani,       # plot the effect of LEN_PORmPUM_LOG:POR_ANIMACY
   type="response", # but plot predicted probabilities
   ylim=c(0, 1),    # ylim: the whole range of possible probabilities
   grid=TRUE)       # & w/ a grid