About the Data

An online survey conducted by Smithson and Shou in 2017 included an experimental component, whereby participants were randomly assigned to receive one of these two questions: “In the next 5,000 years, what is the probability that the human species will become extinct?” or “Within how many years’ time do you expect that the human species will become extinct?”. They also were randomly assigned to make this estimate before or after they had been asked to rank ten existential threats to humanity according to their severity.

The response formats asked respondents to choose a range (e.g., “between 1 in 100 and 1 in 500 chances” vs “100 to 500 years from now”) and then a specific value from a list within that range (e.g., 1 in 300 chances vs 300 years). These alternatives are equated with one another by taking expectations of the probability to arrive at the expected number of years before an extinction event. For instance, an event whose probability is 1 in 300 chances within a given year would be expected to occur once every 300 years.

Both response formats had boundaries on the probabilities that were above 0 and below 1. The lowest probability that could be assigned was less than 1 in 55,000 chances vs more than 55,000 years (.000018), and the highest was more than 1 in 10 chances vs less than 10 years (.1). Thus, these are censored scores.

Samples were taken from three adult populations: The USA (\(N\) = 330), UK (\(N\) = 420), and India (\(N\) = 420). The example here therefore has two experimental factors to include in the model (question format and order of presentation) and one non-experimental factor (nationality).

Model Fitting

This example includes five models, all of them using the logit-logistic distribution with censoring at 0.00001818 and 0.1. The first four models investigate the effects of the two experimental factors (order and format), with the conclusion that order and format have effects in the location submodel but only format does so in the dispersion submodel. The fifth model (m4) includes nationality by using the India sample as the base group, and the coefficients indicate that the US sample rates the probability of human extinction lowest on average, followed by the UK sample.

# Examine effects of format and order for general extinction:
m0 <- cdfquantregC(EQ1_P ~ 1 | 1, fd ='logit',sd ='logistic', c1 = 0.00001818, c2 = 0.1,data = ExtEvent)
m1 <- update(m0, .~. + order + format)
anova(m0, m1)

## Likelihood ratio tests 
## 
##   Resid. Df -2Loglik Df LR stat  Pr(>Chi)    
## 1      1168   -10880                         
## 2      1166   -10912  2  31.741 1.281e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

m2 <- update(m1, .~. |.+ order)
anova(m1, m2)

## Likelihood ratio tests 
## 
##   Resid. Df -2Loglik Df LR stat Pr(>Chi)
## 1      1166   -10912                    
## 2      1165   -10912  1 0.49173   0.4832

m3 <- update(m1, .~. |.+ format)
anova(m1, m3)

## Likelihood ratio tests 
## 
##   Resid. Df -2Loglik Df LR stat Pr(>Chi)   
## 1      1166   -10912                       
## 2      1165   -10921  1   9.169 0.002461 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

m4 <- update(m3, .~. + nation|.)
anova(m3, m4)

## Likelihood ratio tests 
## 
##   Resid. Df -2Loglik Df LR stat  Pr(>Chi)    
## 1      1165   -10921                         
## 2      1163   -10971  2  50.694 9.818e-12 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

summary(m4)

## Family:  logit logistic 
## Call:  cdfquantregC(formula = EQ1_P ~ order + format + nation | format,  
##     data = ExtEvent, fd = "logit", sd = "logistic", c1 = 1.818e-05,  
##     c2 = 0.1) 
## 
## Mu coefficients (Location submodel)
##                Estimate Std. Error z value Pr(>|z|)    
## (Intercept)     -5.7637     0.1548 -37.226  < 2e-16 ***
## orderranklater  -0.4515     0.1411  -3.199 0.001380 ** 
## formatyear      -0.6754     0.1421  -4.752 2.01e-06 ***
## nationUK        -1.1840     0.1651  -7.173 7.33e-13 ***
## nationUS        -0.6288     0.1779  -3.535 0.000408 ***
## 
## Sigma coefficients (Dispersion submodel)
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  0.39009    0.03482  11.204  < 2e-16 ***
## formatyear  -0.14132    0.04894  -2.887  0.00388 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Converge:  successful completion
## Log-Likelihood:  5485.719 
## 
## Gradient:  0

Model Evaluation and Diagnosis

An examination of the parameter estimate correlation matrix reveals no correlations whose magnitudes are alarming. The model distribution is fairly similar to the empirical distribution and the residuals are reasonably well-behaved.

cov2cor(vcov(m4))

##                     (mu)_(Intercept) (mu)_orderranklater (mu)_formatyear
## (mu)_(Intercept)          1.00000000         -0.39421528     -0.54905657
## (mu)_orderranklater      -0.39421528          1.00000000     -0.02406778
## (mu)_formatyear          -0.54905657         -0.02406778      1.00000000
## (mu)_nationUK            -0.49983137         -0.08422835      0.06185522
## (mu)_nationUS            -0.48465439         -0.04338752      0.06946513
## (sigma)_(Intercept)      -0.01524312         -0.02649244      0.01447465
## (sigma)_formatyear        0.01676837          0.03727001     -0.03039308
##                     (mu)_nationUK (mu)_nationUS (sigma)_(Intercept)
## (mu)_(Intercept)      -0.49983137   -0.48465439         -0.01524312
## (mu)_orderranklater   -0.08422835   -0.04338752         -0.02649244
## (mu)_formatyear        0.06185522    0.06946513          0.01447465
## (mu)_nationUK          1.00000000    0.44313215          0.02569131
## (mu)_nationUS          0.44313215    1.00000000          0.01750136
## (sigma)_(Intercept)    0.02569131    0.01750136          1.00000000
## (sigma)_formatyear    -0.03505728   -0.03514061         -0.71236420
##                     (sigma)_formatyear
## (mu)_(Intercept)            0.01676837
## (mu)_orderranklater         0.03727001
## (mu)_formatyear            -0.03039308
## (mu)_nationUK              -0.03505728
## (mu)_nationUS              -0.03514061
## (sigma)_(Intercept)        -0.71236420
## (sigma)_formatyear          1.00000000

plot(m4, nbins = 50)

About the Data

The USA sample in the online survey conducted by Smithson and Shou (2017) as described earlier included items from the social and economic conservatism scales created by Everett (2013). Each item asked respondents to rate their feelings about the issue described in the item on a scale from 0 to 100, according to this instruction: “Please indicate the extent to which you feel positive or negative towards each issue. Scores of 0 indicate greater negativity, and scores of 100 indicate greater positivity. Scores of 50 indicate that you feel neutral about the issue.”

Model Fitting

The next figure shows a histogram of the ratings on the issue of “gun ownership”. This is clearly a strongly polarizing issue. There are reasonable arguments for treating the bounds on the gun ownership scale either as censored scores or true scores. Here, we treat the bounds as true scores, so that responses are considered as a doubly-bounded random variable.

Histograms of gun ownership ratings separated by political orientation show clear differences among the four orientations. The sources of the polarization in the distribution are primarily the Democrats and Republicans, as would be expected. We should expect an accurate model to highlight these differences, given that there are sufficiently many people in each of the four groups for such a model to detect sizable group differences.

# How many people occupy the political orientation groups in the sample?
table(gunowndata$political)

## 
##    Democrat Independent      NoPref  Republican 
##          98         109          46          68

# 
par(mfrow = c(2,2),mar = c(4,4,1,1))
truehist(gunowndata$gunown[gunowndata$political == "Democrat"], nbins = 50, main = "Democrat", xlab = "gun ownership", ylab = "density", ylim = c(0,11), col = "red")
truehist(gunowndata$gunown[gunowndata$political == "Independent"], nbins = 50, main = "Independent", xlab = "gun ownership", ylab = "density", ylim = c(0,11), col = "red")
truehist(gunowndata$gunown[gunowndata$political == "NoPref"], nbins = 50, main = "No Preference", xlab = "gun ownership", ylab = "density", ylim = c(0,11), col = "red")
truehist(gunowndata$gunown[gunowndata$political == "Republican"], nbins = 50, main = "Republican", xlab = "gun ownership", ylab = "density", ylim = c(0,11), col = "red")

The first three models test for the effect of political orientation in the non-hurdle component of the data, using the burr8-burr8 distribution. Including political orientation in the dispersion submodel does not improve model fit, so subsequent models omit it.

mod0 <- cdfquantregH(gunown ~ 1, zero.fo = ~1, one.fo = ~1, fd = 'burr8', sd = 'burr8', type = 'ZO', data = gunowndata)
mod1 <- cdfquantregH(gunown ~ political, zero.fo = ~1, one.fo = ~1, fd = 'burr8', sd = 'burr8', type = 'ZO', data = gunowndata)
mod2 <- cdfquantregH(gunown ~ political|political, zero.fo = ~1, one.fo = ~1, fd = 'burr8', sd = 'burr8', type = 'ZO', data = gunowndata)
mod3 <- cdfquantregH(gunown ~ political, zero.fo = ~political, one.fo = ~political, fd = 'burr8', sd = 'burr8', type = 'ZO', data = gunowndata)


anova(mod1,mod3)

## Likelihood ratio tests 
## 
##   Resid. Df -2Loglik Df LR stat Pr(>Chi)    
## 1       314   353.83                        
## 2       308   322.05  6  31.777  1.8e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

summary(mod3)

## Family:  burr8 burr8 
## Call:  cdfquantregH(formula = gunown ~ political, data = gunowndata,  
##     fd = "burr8", sd = "burr8", zero.fo = ~political, one.fo = ~political,  
##     type = "ZO") 
## 
## Mu coefficients (Location submodel)
##                      Estimate Std. Error z value Pr(>|z|)    
## (Intercept)           -0.2785     0.1296  -2.148 0.031677 *  
## politicalIndependent   0.6523     0.1796   3.633 0.000281 ***
## politicalNoPref        0.2917     0.2340   1.247 0.212569    
## politicalRepublican    1.0886     0.2020   5.390 7.03e-08 ***
## 
## Sigma coefficients (Dispersion submodel)
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -0.19454    0.05573  -3.491 0.000482 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero component coefficients
##                      Estimate Std. Error z value Pr(>|z|)    
## (Intercept)           -1.7918     0.2887  -6.207 5.41e-10 ***
## politicalIndependent  -1.4759     0.5855  -2.521   0.0117 *  
## politicalNoPref        0.5108     0.4595   1.112   0.2662    
## politicalRepublican   -1.7047     0.7736  -2.204   0.0276 *  
## 
## One component coefficients
##                      Estimate Std. Error z value Pr(>|z|)    
## (Intercept)           -3.8712     0.7144  -5.419    6e-08 ***
## politicalIndependent   1.6841     0.7820   2.154  0.03126 *  
## politicalNoPref        1.5198     0.8855   1.716  0.08611 .  
## politicalRepublican    2.3308     0.7820   2.980  0.00288 ** 
## 
## Converge:  
## Log-Likelihood:  -161.0249

The final model shows the expected effects of political orientation in all three model components. The location submodel yields higher ratings for Republicans and Independents than for Democrats, whereas the submodel does not find a significant difference between the Democrat and No Preference groups. These differences are echoed in the zero and one components. Republicans and Independents are more likely to give zero ratings and less likely to give ratings of one than Democrats. The No Preference group has a marginally greater tendency than Democrats to give ratings of 1, but it does not reach significance.