表2.2の再現

Ross2012のデータを読み込む。

load("Replication data for The Oil Curse - Ross 2012.RData")

浜中によるコードは以下の通り。

表2.1と同じく、浜中コードでは非産油国と産油国の順番が入れ替わっている点に注意が必要。

library(dplyr)
TAXa <- filter(x, year==2002)
TAXaa <- filter(TAXa, gdpcap2000_sup<5000)
TAXab <- filter(TAXa, gdpcap2000_sup>5000)
t.test(TAXaa$tax_goods[TAXaa$oil_gas_value100==1],
       TAXaa$tax_goods[TAXaa$oil_gas_value100==0])
t.test(TAXab$tax_goods[TAXab$oil_gas_value100==1],
       TAXab$tax_goods[TAXab$oil_gas_value100==0])
t.test(TAXa$tax_goods[TAXa$oil_gas_value100==1],
       TAXa$tax_goods[TAXa$oil_gas_value100==0])

1行目(低所得国)のt検定


    Welch Two Sample t-test

data:  TAXaa$tax_goods[TAXaa$oil_gas_value100 == 1] and TAXaa$tax_goods[TAXaa$oil_gas_value100 == 0]
t = -2.0134, df = 16.947, p-value = 0.06024
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -16.1385649   0.3793741
sample estimates:
mean of x mean of y 
 24.91274  32.79234 

平均値は表と同じだが、有意水準が異なる。浜中コードでは有意水準は10%だが(p値が0.06)、表では5%水準となっており、ロスの分析ではp値は0.0371である。また表2.2のt検定は片側だが、浜中コードは両側検定になっているので、ここを修正する必要がある。また、ロスは等分散を前提としているので、これを再現実験に反映させる必要がある。さらに表2.1で行なったように、“t.test”の書式をシンプルなものに変更する。これらの点を考慮すると、コードは以下のように修正できる。

t.test(tax_goods~oil_gas_value100,data=TAXaa,
       alternative="greater",var.equal = T)

    Two Sample t-test

data:  tax_goods by oil_gas_value100
t = 1.8193, df = 57, p-value = 0.03706
alternative hypothesis: true difference in means is greater than 0
95 percent confidence interval:
 0.637876      Inf
sample estimates:
mean in group 0 mean in group 1 
       32.79234        24.91274 

これでp値は.0371となり、ロスの分析結果と同じになった。

2行目(低所得国)のt検定


    Welch Two Sample t-test

data:  TAXab$tax_goods[TAXab$oil_gas_value100 == 1] and TAXab$tax_goods[TAXab$oil_gas_value100 == 0]
t = -1.3831, df = 14.784, p-value = 0.1872
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -14.64076   3.12651
sample estimates:
mean of x mean of y 
 24.08763  29.84475 

浜中コードの結果は、表とは異なる。平均値が若干異なり、p値も異なる。ロスの分析結果では、p値は0.070である。ロスは等分散を前提に分析を行なっており、また片側検定なので、これに合わせてコードを以下のように修正する(“t.test”の書式も変更済み)。

t.test(tax_goods~oil_gas_value100,data = TAXab,
       alternative = "greater",var.equal = T)

    Two Sample t-test

data:  tax_goods by oil_gas_value100
t = 1.5681, df = 35, p-value = 0.06293
alternative hypothesis: true difference in means is greater than 0
95 percent confidence interval:
 -0.4459468        Inf
sample estimates:
mean in group 0 mean in group 1 
       29.84475        24.08763 

ロスの分析結果に近づくが、これでもやや異なる。原因は非産油国の平均値が異なる点にあり、ロスの分析では29.55708だが、浜中コードと上記修正コードでは29.84475である。これはサンプルの取り扱いが異なっていることを示しており、表2.1と同じ問題が発生していることをうかがわせる。

Rでサンプル数とNA数を確認する。

table(is.na(TAXab$tax_goods),TAXab$oil_gas_value100)
       
         0  1
  FALSE 26 11
  TRUE   3  7

非産油国のNAは3(計29)、産油国のNAは7(計18)、NAは計10、合計47である。ロスがSTATAで行なった分析では、非産油国が27、産油国が11であり、やはりサンプル数が異なる。この点については表2.1と同様、STATAの書類を確認しないといけないので、現段階ではこのままとする。

3行目(すべて)のt検定


    Welch Two Sample t-test

data:  TAXa$tax_goods[TAXa$oil_gas_value100 == 1] and TAXa$tax_goods[TAXa$oil_gas_value100 == 0]
t = -2.5118, df = 35.167, p-value = 0.01675
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -12.887169  -1.367755
sample estimates:
mean of x mean of y 
 24.50018  31.62765 

浜中コードでは5%水準で有意(p値が0.018)だが、表では1%水準で有意となっている(ロスの分析ではp値は0.0077)。ロスの分析に合わせて、片側検定で等分散を前提とすると、浜中コードはロスの分析結果に一致する(“t.test”の書式も変更済み)。

t.test(tax_goods~oil_gas_value100,data = TAXa,
       alternative = "greater",var.equal = T)

    Two Sample t-test

data:  tax_goods by oil_gas_value100
t = 2.4698, df = 95, p-value = 0.007652
alternative hypothesis: true difference in means is greater than 0
95 percent confidence interval:
 2.333925      Inf
sample estimates:
mean in group 0 mean in group 1 
       31.62765        24.50018 

これでロスの分析結果と同じになる。