Did Wongupparaj et al. debunk the Flynn Effect?

Wongupparaj et al. (2017) released a study that found that backward digit span/corsi span scores declined, while forward digit span/corsi scores increased. The idea is that backward memory is more g-loaded than forward memory, so this supports the view that the Flynn Effect is an artefact and that true IQ is lowering.

This wrong. Let me quickly pull up their dataset.

head(savy)

Recreating the first regression:

qwe <- lm(data=savy, Mean ~ Year, weights=Size)
summary(qwe)

Call:
lm(formula = Mean ~ Year, data = savy, weights = Size)

Weighted Residuals:
    Min      1Q  Median      3Q     Max 
-58.915  -4.623   0.428   4.564  51.829 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)  
(Intercept) 16.109773   8.273144   1.947    0.052 .
Year        -0.005778   0.004121  -1.402    0.161  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 9.877 on 592 degrees of freedom
Multiple R-squared:  0.003309,  Adjusted R-squared:  0.001626 
F-statistic: 1.966 on 1 and 592 DF,  p-value: 0.1614

Alright. Replicated. What about their second regression?

qwe <- lm(data=savy, Mean ~ Year + Age + Male + Type_C + Platform, weights=Size)
summary(qwe)

Call:
lm(formula = Mean ~ Year + Age + Male + Type_C + Platform, data = savy, 
    weights = Size)

Weighted Residuals:
    Min      1Q  Median      3Q     Max 
-49.643  -4.525   0.440   4.506  51.708 

Coefficients:
              Estimate Std. Error t value        Pr(>|t|)    
(Intercept)  7.7411426  8.0485412   0.962           0.337    
Year        -0.0020717  0.0039988  -0.518           0.605    
Age          0.0058195  0.0014088   4.131 0.0000414092464 ***
Male        -0.0008939  0.0021264  -0.420           0.674    
Type_C       0.8781821  0.1312821   6.689 0.0000000000523 ***
Platform    -0.0862458  0.1244160  -0.693           0.488    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 9.418 on 588 degrees of freedom
Multiple R-squared:  0.09989,   Adjusted R-squared:  0.09223 
F-statistic: 13.05 on 5 and 588 DF,  p-value: 0.000000000004612

Different from their 2A table. Not sure what’s going on. Year is not statistically signfiicant anymore.

---

So, that aside, I labelled all 598 effect sizes by country of origin (22 cases could not be identified) and socioeconomic status (factor with levels 0, 1, 2) to test for omitted variable biases. I also corrected a few clerical errors, and added a few studies (6?) I randomly found when searching for the publications for each effect size.

subset <- testtempyarfwefe %>% filter(type=='backward' & ((!Country==1 & !Country==0 & !Country=='') | is.na(Country)))
subset$Country[is.na(subset$Country)] <- 'Unknown'
subset$SES[is.na(subset$SES)] <- 1
subset$SE <- subset$SD/sqrt(subset$Size)

Then, I noticed that the mean and the standard deviation of backward digit span was non-linearly related to age, even after controlling for specific country and SES, so I normed scores using splines instead of a linear method.

p <- splinechart(datos=subset, var1='Age', var2='Mean', knots=5, breakn=10, metastyle=F, label='Age', labely='Mean Backward Digit Span',
                 covariatefunc=' + as.factor(SES) + Country', ws='Size', burndata=data.frame(SES=1, Country='UK'), outputyes=F)
p

#Yes, this does pass significance testing
p <- splinechart(datos=subset, var1='Age', var2='SD', knots=5, breakn=10, metastyle=F, label='Age', labely='SD Backward Digit Span',
                 covariatefunc=' + as.factor(SES) + Country', ws='Size', burndata=data.frame(SES=1, Country='UK'), outputyes=F)
p

Onto the hard calculations. Now I have scores normed for age, both in terms of variance and mean.

lr <- lm(data=subset, Mean ~ rcs(Age, 5) + as.factor(SES), weights=sqrt(Size))
summary(lr)

Call:
lm(formula = Mean ~ rcs(Age, 5) + as.factor(SES), data = subset, 
    weights = sqrt(Size))

Weighted Residuals:
     Min       1Q   Median       3Q      Max 
-10.2619  -1.2523  -0.2109   1.0086   7.9391 

Coefficients:
                  Estimate Std. Error t value             Pr(>|t|)    
(Intercept)        0.59841    0.24813   2.412             0.016181 *  
rcs(Age, 5)Age     0.31503    0.02041  15.436 < 0.0000000000000002 ***
rcs(Age, 5)Age'   -5.34306    0.47150 -11.332 < 0.0000000000000002 ***
rcs(Age, 5)Age''   7.45751    0.69468  10.735 < 0.0000000000000002 ***
rcs(Age, 5)Age''' -2.20123    0.26302  -8.369 0.000000000000000419 ***
as.factor(SES)1    0.54224    0.15892   3.412             0.000689 ***
as.factor(SES)2    0.52464    0.18372   2.856             0.004446 ** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.214 on 592 degrees of freedom
Multiple R-squared:  0.4285,    Adjusted R-squared:  0.4228 
F-statistic: 73.99 on 6 and 592 DF,  p-value: < 0.00000000000000022

metaobj2 <- metafor::rma(yi=Mean, sei=SE, data=subset, mods = ~ rcs(Age, 5) + as.factor(SES))
summary(metaobj2)

Mixed-Effects Model (k = 599; tau^2 estimator: REML)

   logLik   deviance        AIC        BIC       AICc   
-681.8482  1363.6963  1379.6963  1414.7644  1379.9433   

tau^2 (estimated amount of residual heterogeneity):     0.5408 (SE = 0.0339)
tau (square root of estimated tau^2 value):             0.7354
I^2 (residual heterogeneity / unaccounted variability): 98.06%
H^2 (unaccounted variability / sampling variability):   51.42
R^2 (amount of heterogeneity accounted for):            48.41%

Test for Residual Heterogeneity:
QE(df = 592) = 35716.9755, p-val < .0001

Test of Moderators (coefficients 2:7):
QM(df = 6) = 525.5336, p-val < .0001

Model Results:

---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

lr$coefficients <- metaobj2$b

metaobj32 <- metafor::rma(yi=SD, sei=SE, data=subset, mods = ~ rcs(Age, 5) + as.factor(SES))
summary(metaobj2)

Mixed-Effects Model (k = 599; tau^2 estimator: REML)

   logLik   deviance        AIC        BIC       AICc   
-681.8482  1363.6963  1379.6963  1414.7644  1379.9433   

tau^2 (estimated amount of residual heterogeneity):     0.5408 (SE = 0.0339)
tau (square root of estimated tau^2 value):             0.7354
I^2 (residual heterogeneity / unaccounted variability): 98.06%
H^2 (unaccounted variability / sampling variability):   51.42
R^2 (amount of heterogeneity accounted for):            48.41%

Test for Residual Heterogeneity:
QE(df = 592) = 35716.9755, p-val < .0001

Test of Moderators (coefficients 2:7):
QM(df = 6) = 525.5336, p-val < .0001

Model Results:

---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

lr3 <- lm(data=subset, SD ~ rcs(Age, 5) + as.factor(SES), weights=sqrt(Size))
summary(lr)

Call:
lm(formula = Mean ~ rcs(Age, 5) + as.factor(SES), data = subset, 
    weights = sqrt(Size))

Weighted Residuals:
     Min       1Q   Median       3Q      Max 
-10.2619  -1.2523  -0.2109   1.0086   7.9391 

Coefficients:
     Estimate Std. Error t value             Pr(>|t|)    
[1,]  0.57151    0.24813   2.303              0.02161 *  
[2,]  0.30312    0.02041  14.852 < 0.0000000000000002 ***
[3,] -4.99721    0.47150 -10.599 < 0.0000000000000002 ***
[4,]  6.95332    0.69468  10.009 < 0.0000000000000002 ***
[5,] -2.03386    0.26302  -7.733   0.0000000000000455 ***
[6,]  0.62281    0.15892   3.919   0.0000992599476217 ***
[7,]  0.57459    0.18372   3.128              0.00185 ** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.214 on 592 degrees of freedom
Multiple R-squared:  0.4285,    Adjusted R-squared:  0.4228 
F-statistic: 73.99 on 6 and 592 DF,  p-value: < 0.00000000000000022

lr3$coefficients <- metaobj32$b

uzi <- seq(from=4, to=92.5, by=.1)
uzi2 <- data.frame(Age=uzi, SES=1)
uzi2$SD_age = predict(lr3, uzi2)

subset$residualscore <- lr$residuals
subset$roundage <- format(round(subset$Age, 1), nsmall = 1)
uzi2$Age_formatted <- format(uzi2$Age, nsmall = 1)

subset2 <- left_join(
  subset,
  uzi2 %>% select(Age_formatted, SD_age),
  by = c('roundage' = 'Age_formatted')
)

subset2$IQ <- subset2$residualscore/subset2$SD_age*15+100

Do these scores normed for age/SES have a Flynn Effect? Yes. See the regression table – increase of 0.14 (SE = 0.05) points every year.

p <- splinechart(datos=subset2, var1='Year', var2='IQ', knots=5, breakn=5, metastyle=F, label='Year', labely='IQ',
                 covariatefunc='+ Country', ws='Size', burndata=data.frame(SES=1, Country='UK'), outputyes = T)

Call:
lm(formula = paste0("var22", modell), data = datos, weights = sqrt(defaultweights))

Weighted Residuals:
    Min      1Q  Median      3Q     Max 
-133.79  -14.29   -1.38   11.05  121.40 

Coefficients:
                      Estimate Std. Error t value Pr(>|t|)   
(Intercept)         -183.98760   98.45329  -1.869  0.06218 . 
var11                  0.13720    0.04884   2.809  0.00514 **
CountryAustralia      15.13592    5.55526   2.725  0.00664 **
CountryAustria        27.19053   13.96176   1.948  0.05197 . 
CountryBelgium         5.27791    5.73316   0.921  0.35766   
CountryBrazil          3.67153    5.89434   0.623  0.53361   
CountryCanada         15.09410    5.61930   2.686  0.00744 **
CountryChina          11.66256    6.24306   1.868  0.06227 . 
CountryColombia       13.18377   12.80726   1.029  0.30374   
CountryCroatia         5.77719    7.12118   0.811  0.41756   
CountryDenmark        -1.42541    8.03405  -0.177  0.85924   
CountryFinland        12.67420    6.06690   2.089  0.03715 * 
CountryFrance         14.56088    6.94555   2.096  0.03649 * 
CountryGermany         8.75160    5.44791   1.606  0.10875   
CountryHong Kong      17.75716    7.34317   2.418  0.01592 * 
CountryHungary        14.17243   10.14527   1.397  0.16298   
CountryIndia          -7.70436    8.91570  -0.864  0.38788   
CountryIran            1.48125    9.76069   0.152  0.87943   
CountryIsrael          8.57904    7.87656   1.089  0.27654   
CountryItaly           3.84550    5.33187   0.721  0.47107   
CountryJamaica         2.43182    9.47772   0.257  0.79759   
CountryJapan           9.85085    6.90302   1.427  0.15413   
CountryMexico         -4.90161    5.95809  -0.823  0.41104   
CountryNetherlands    10.68805    5.54619   1.927  0.05447 . 
CountryNew Zealand     2.77731    5.75294   0.483  0.62945   
CountryNorway          8.04522    6.86588   1.172  0.24179   
CountryPortugal        2.14251    6.50921   0.329  0.74217   
CountrySingapore       4.85551    6.45849   0.752  0.45249   
CountrySlovenia       24.07165   11.62842   2.070  0.03890 * 
CountrySouth Africa   -9.97060    8.28712  -1.203  0.22943   
CountrySouth Korea     2.88260    6.19272   0.465  0.64177   
CountrySpain          11.49762    5.57165   2.064  0.03952 * 
CountrySweden          7.43187    5.53835   1.342  0.18017   
CountryTaiwan         -1.40316    6.98551  -0.201  0.84088   
CountryTunisia        10.32997    6.80912   1.517  0.12981   
CountryTurkey          9.78593    9.74573   1.004  0.31575   
CountryUK              8.98739    5.51222   1.630  0.10357   
CountryUnknown        10.78941    5.41159   1.994  0.04666 * 
CountryUSA            10.15682    5.20103   1.953  0.05134 . 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 27.44 on 560 degrees of freedom
Multiple R-squared:  0.1803,    Adjusted R-squared:  0.1247 
F-statistic: 3.241 on 38 and 560 DF,  p-value: 0.000000001205

[1] "Non-linear vs linear:"
Analysis of Variance Table

Model 1: var22 ~ rcs(var11, 5) + Country
Model 2: var22 ~ var11 + Country
  Res.Df    RSS Df Sum of Sq      F Pr(>F)
1    557 420168                           
2    560 421776 -3   -1607.9 0.7105  0.546

p

Same attempt, but with cohort effects instead. No increase or decrease.

subset2$cohort <- subset2$Year - subset2$Age

p <- splinechart(datos=subset2, var1='cohort', var2='IQ', knots=5, breakn=10, metastyle=F, label='Cohort', labely='IQ',
                 covariatefunc='+ Country', ws='Size', burndata=data.frame(SES=1, Country='UK'), outputyes = F)
p

Calculation of international differences (resetting data, controlling for year). Obviously, these means don’t have a reference population, and I’m sure that a dozen criticisms could be hauled at them, but I figured people would be curious anyway, so why not plot them.

subset <- testtempyarfwefe %>% filter(type=='backward' & ((!Country==1 & !Country==0 & !Country=='') | is.na(Country)))
subset$Country[is.na(subset$Country)] <- 'Unknown'
subset$SES[is.na(subset$SES)] <- 1
subset$SE <- subset$SD/sqrt(subset$Size)


########
lr <- lm(data=subset, Mean ~ rcs(Age, 5) + as.factor(SES) + Year, weights=sqrt(Size))
metaobj2 <- metafor::rma(yi=Mean, sei=SE, data=subset, mods = ~ rcs(Age, 5) + as.factor(SES) + Year)
lr$coefficients <- metaobj2$b

metaobj32 <- metafor::rma(yi=SD, sei=SE, data=subset, mods = ~ rcs(Age, 5) + as.factor(SES) + Year)
lr3 <- lm(data=subset, SD ~ rcs(Age, 5) + as.factor(SES) + Year, weights=sqrt(Size))
lr3$coefficients <- metaobj32$b

uzi <- seq(from=4, to=92.5, by=.1)
uzi2 <- data.frame(Age=uzi, SES=1, Year=2008)
uzi2$SD_age = predict(lr3, uzi2)

subset$residualscore <- lr$residuals
subset$roundage <- format(round(subset$Age, 1), nsmall = 1)
uzi2$Age_formatted <- format(uzi2$Age, nsmall = 1)

subset2 <- left_join(
  subset,
  uzi2 %>% select(Age_formatted, SD_age),
  by = c('roundage' = 'Age_formatted')
)

subset2$IQ <- subset2$residualscore/subset2$SD_age*15+100
subset2$pasteID <- paste0(subset2$Authors, subset2$Year)

zxc <- unique(subset2 %>% group_by(pasteID) %>% summarise(IQ2 = mean(IQ), n = sum(Size), Country = Country))

Warning: Returning more (or less) than 1 row per `summarise()` group was deprecated in dplyr 1.1.0.
Please use `reframe()` instead.
When switching from `summarise()` to `reframe()`, remember that `reframe()` always returns an ungrouped data frame and adjust accordingly.`summarise()` has grouped output by 'pasteID'. You can override using the `.groups` argument.

zxc$SE <- 1/sqrt(zxc$n)

zxcv <- unique(zxc %>% group_by(Country) %>% summarise(unwtdIQ = mean(IQ2), wtdIQ = wtd.mean(IQ2, n), medianIQ = median(IQ2), n2 = sum(n), k = n()))
zxcv  %>% arrange(-medianIQ)