Niet beter dan de Heckman two-step : een Monte-Carlo-onderzoek naar de Sample-selectionmethode van Newey

(1)

Niet beter dan de Heckman two-step

Een Monte-Carlo-onderzoek naar de Sample-Selectionmethode van Newey

Bernd Schouten Studentnummer: 10356851

Bachelorscriptie onder begeleiding van Dr. J.C.M. van Ophem

Faculteit Economie en Bedrijfskunde Econometrie en Operationele Research

(2)

Verklaring eigen werk

Hierbij verklaar ik, Bernd Schouten, dat ik deze scriptie zelf geschreven heb en dat ik de vol-ledige verantwoordelijkheid op me neem voor de inhoud ervan. Ik bevestig dat de tekst en het werk dat in deze scriptie gepresenteerd wordt origineel is en dat ik geen gebruik heb gemaakt van andere bronnen dan die welke in de tekst en in de referenties worden genoemd. De Faculteit Economie en Bedrijfskunde is alleen verantwoordelijk voor de begeleiding tot het inleveren van de scriptie, niet voor de inhoud.

Excerpt In this study, I investigate the differences between two methods for estimating samples with sample selection bias: the Heckman two-step method and the alternative proposed by Newey (2009). Through computer simulation I study what happens when the correlation between the error terms and the correlation between the independent variables. The results from the Monte Carlo simulation indicate that the differences between the results of the methods are small.

1 Inleiding

Veel econometrische datasets hebben last van sample selection bias, wat wil zeggen dat de ge-schatte waarden een vertekend beeld geven. Een gebruikelijke manier om dit uit te leggen is met het volgende, zeer beroemde, voorbeeld uit de economie (zo bijvoorbeeld Puhani, 2000): men wil de hoogte van het loon van mensen verklaren met de hoeveelheid scholing die ze geno-ten hebben. Een relatief groot aantal mensen met weinig scholing werkt echter helemaal niet, omdat het hun aangeboden loon lager is dan hun reserveringsloon. Dat heeft als gevolg dat de mensen met weinig scholing die wel werken een relatief hoog aangeboden loon hebben gekre-gen. Een gewone kleinstekwadratenschatter zal dus de negatieve effecten van weinig scholing onderschatten. Daarin zit het probleem van de sample selection bias.

Dit probleem komt in veel economische onderzoeken voor, ook in minder theoretische situ-aties. Er is bijvoorbeeld een interessant onderzoek gedaan naar de effecten van subsidies van de Duitse overheid op het gebied van Research en Development (Hussinger, 2008). De over-heid geeft namelijk subsidies aan bedrijven om innovatie te stimuleren, maar het is natuurlijk goed mogelijk dat bedrijven zonder de subsidie ook het onderzoek zouden doen, maar het dan zelf zouden financieren. Het is echter niet zinvol om alleen te kijken of bedrijven die subsidie krijgen meer geld uitgeven aan innovatie dan bedrijven die de subsidie niet krijgen. De Duitse

(3)

overheid deelt niet willekeurig geld uit, maar geeft haar subsidiegeld alleen aan veelbelovende bedrijven. Ook hier is dus naar alle waarschijnlijkheid sprake van sample selection bias.

Om toch bruikbare schattingsresultaten te krijgen, is al in 1976 een oplossing bedacht: de Heckman tweestapsmethode (Heckman, 1976, 1979), die ook de Limited Information Maxi-mum Likelihood (LIML) genoemd wordt. Deze methode is relatief eenvoudig en wordt daarom veel gebruikt in empirisch onderzoek, maar dat is niet zo vanzelfsprekend als het lijkt. Heckman zelf (1979, p. 160) merkt namelijk al op dat zijn methode vooral gebruikt zou moeten worden als eerste indruk en niet als definitieve conclusie. Er is dan ook al veel onderzoek gedaan om betere methodes te vinden en die te vergelijken met de LIML.

Een van de alternatieven voor Heckman is de methode van Newey (2009), die enkele van de grootste problemen van de LIML oplost. Als deze methode werkt zoals verwacht wordt, is ze namelijk toepasbaar onder veel minder strenge voorwaarden dan Heckmans schatter. In dit onderzoek wordt de methode van Newey vergeleken met die van Heckman onder verschillende omstandigheden, om te kijken welke van de twee de beste schattingsresultaten geeft.

Daarvoor wordt eerst de Heckman-methode uitgelegd en wordt een aantal van de grootste bezwaren ervan besproken. Vervolgens wordt het alternatief van Newey uitgelegd. Op basis daarvan wordt een Monte-Carlosimulatie ontworpen en uitgevoerd waarin deze twee methodes met elkaar vergeleken worden.

2 Theoretische achtergrond

2.1 De methode van Heckman

De tweestapsmethode (Heckman, 1979) wordt in dit hoofdstuk uitgelegd in haar eenvoudigste vorm. Het is natuurlijk mogelijk om dit te generaliseren naar meer selectievergelijkingen, maar het principe blijft hetzelfde, dus is het handiger om het model zo klein mogelijk te houden. Hier-bij volg ik de notatie en volgorde van uitleggen van Puhani (2000), omdat die overzichtelijker is dan die van Heckman zelf.

(4)

Het model dat geschat moet worden, ziet er als volgt uit voor elke ide waarneming:

y∗_1i = x_1i0 β1+ u1i (1)

y∗_2i = x_2i0 β2+ u2i (2)

y1i = y∗1i als y∗2i> 0 (3)

y1i = 0 als y∗2i≤ 0 (4)

Er zijn goede theoretische gronden om er van uit te gaan dat de storingstermen u1i en u2i niet

ongecorreleerd zijn. In het voorbeeld van werkloosheid uit de inleiding vertaalt zich dat als volgt: mensen met weinig scholing die vanwege niet waarneembare redenen verwachten dat ze weinig geld zullen verdienen (dus een relatief lage u1i hebben), zullen minder gemotiveerd

zijn om een baan te zoeken (en hebben dus ook een lage u2i). Er moet bij het schatten dus

rekening gehouden worden met het feit dat de storingen positief gecorreleerd zijn. Wat verder nog aangenomen wordt, is dat de storingen bivariaat normaal verdeeld zijn. Dat is natuurlijk een strenge aanname en een van de grootste bezwaren tegen deze methode. Het zorgt ervoor dat de storingen er als volgt uitzien:

  u₁ u2  ∼ BN     0 0  ,   σ₁2 σ12 σ12 σ₂2     (5)

Met deze modelaannamen is het mogelijk om een maximum-likelihoodschatter op te stellen. Het berekenen van deze Full Information Maximum Likelihood (FIML) kost echter veel com-puterkracht. Dit was zeker in de jaren zeventig toen Heckman zijn methode voorstelde een groot probleem. Daarom was er dus behoefte aan een eenvoudiger te berekenen alternatief.

Heckman suggereerde dat uiteindelijk deze vergelijking geschat moet worden:

E(y∗_1i|x1i, y∗2i> 0) = x01iβ1+ E(u1i|u2i> −x02iβ2) (6)

Hier is ook goed zichtbaar dat als er geen correlatie is tussen de storingen dit de gewone klein-stekwadratenschatter is. Dan is immers de verwachte waarde E(u1i|u2i> −x02iβ2) gelijk aan de

verwachte waarde van u1, waarbij alle termen x2en u2zijn weggevallen. Onder de aanname van

normaliteit uit vergelijking (5) is het mogelijk om de conditionele verwachting van u1expliciet

(5)

E(u1i|u2i> −x02iβ2) = σ_σ12₂ φ (−x

0 2iβ2/σ2)

1−Φ(−x0_2iβ2/σ2) (7)

Waarbij, zoals gebruikelijk, φ (..) en Φ(..) de standaardnormale dichtheids- en verdelingsfunctie uitdrukken.

De eerste stap van de tweestapsmethode bestaat er uit om de zogeheten inverse Mills-ratio λ te schatten met behulp van een Probit-model (Voor het Probit-model, zie bijvoorbeeld Heij et al. (2004) pp. 444-452). λ (−x0_2iβˆ2/ ˆσ2) = φ (−x 0 2iβˆ2/ ˆσ2) 1−Φ(−x0_2iβˆ2/ ˆσ2) (8)

De tweede stap is dan om met de geschatte waarde voor de inverse Mills-ratio, vergelijking (6) te schatten. Daarbij wordt σ12

σ2 als een van de te schatten parameters toegevoegd en de volgende

regressievergelijking geschat: E(y∗_1i|x1i, y∗2i> 0) = x1i0 β1+ σ12 σ2 λ (−x 0 2i ˆ β2 σ2) + εi (9)

Hierbij is εieen storingsterm met verwachting 0 en variantie V (εi) = σ12− σ₁₂2 σ₂2 h_x0 2iβ σ2 λ ( x0_2iβ σ2 ) + λ ( x0_2iβ σ2 ) 2i

Het is eenvoudig te zien dat εi heteroskedastisch is, aangezien zijn variantie afhankelijk is van

x₂. Dat zorgt ervoor dat de tweestapsschatter niet efficiënt is.

De kracht van de methode zit erin dat het probleem hier als een bijzondere vorm van een omitted variable biaswordt beschouwd, waarbij de inverse Mills-ratio de omitted variable is. Het is dus van een relatief ingewikkelde maximum-likelihoodschatter naar een vrij overzich-telijk te schattten probleem teruggebracht. Dat is de reden dat deze methode zeer veelvuldig gebruikt is. Er kleven echter ook veel nadelen aan de tweestapsmethode.

2.2 Problemen met de Heckman-methode

Het belangrijkste nadeel van de LIML en de FIML is dat de variabelen die de selectie verkla-ren x2 vaak sterk overlappen met de verklarende variabelen x1 of zelfs helemaal identiek zijn

daaraan. Het enige dat er in het laatste geval voor zorgt dat er geen problemen met multicolli-neariteit ontstaan is de non-limulticolli-neariteit van de inverse Mills-ratio. Die ratio is echter in een groot

(6)

deel van zijn domein wel (bijna) lineair, zoals duidelijk te zien is in de onderstaande grafiek:

Figuur 1 De quasi-lineariteit van de inverse Mills-ratio (Puhani, 2000)

De logische oplossing hiervoor is variabelen te vinden in x2 die voldoende zeggen over y2,

maar niet of nauwelijks gecorreleerd zijn met y1, zoals wordt voorgesteld door onder andere

Little en Rubin (1987, p. 230). In de praktijk bestaan zulke variabelen vaak niet. Dan is er dus sprake van grote mate van multicollineariteit en is de schatter daardoor bijzonder instabiel.

Het is op puur theoretische gronden echter niet mogelijk om te bepalen hoe groot de gevol-gen precies zijn. Daarom is er al veel onderzoek gedaan met behulp van Monte-Carlosimulaties. Het grote voordeel van computersimulaties is natuurlijk dat de techniek wordt toegepast op vari-abelen waarvan de verdeling en onderliggende datagenererend proces bekend is. Nawata (1993) richt zich bijvoorbeeld op de verschillen tussen de Heckman tweestapsmethode en OLS. Hij wil onderzoeken wat er gebeurt als de correlatie tussen de afhankelijke variabelen (x1 en

x₂) en tussen de storingen (u1en u2) groter worden. Hij gebruikt daarvoor het volgende model

(7)

x1i = ξ1i (10) x2i = [ρ1ξ1i+ (1 − ρ1)ξ2i]/ q ρ₁2+ (1 − ρ1)2 (11) u1i = ε1i (12) u2i = [ρ210ε1i+ (1 − ρ2)ε2i]/ q ρ₂2+ (1 − ρ2)2 (13)

Hierbij zijn de ξ ’s i.i.d. verdeeld uit een uniforme (0,20] verdeling en de ε’s i.i.d. uit een bivari-aat normale verdeling met verwachting 0 en variantie 1. Wat Nawata vervolgens onderzoekt, is wat er gebeurt als ρ1en ρ2verschillende waarden aannemen. De verwachting is dat als één van

die twee groot is, de LIML-schatter een groter bias krijgt en de variantie ervan groot zal zijn. Het interessante wat hij ontdekt is echter dat dat alleen gebeurt als ρ1, de correlatie tussen x1en

x2, groot is (vooral ≥ 0.9), onafhankelijk van hoe groot ρ2is. Als ρ2groot is en ρ1klein, dan is

de bias van de OLS-schattingen groot, maar doet de LIML het nog prima. Alleen de correlatie tussen de verklarende variabelen is dus van belang voor de prestatie van de Heckman-methode.

2.3 Alternatieve schatters

Een andere reden waarom Heckman minder gebruikt wordt, is omdat er inmiddels andere me-thoden zijn die minder strenge aannames vereisen. De meest voor de hand liggende is de ge-wone, Full Information Maximum Likelihood (FIML). In 1976, toen Heckman zijn methode bedacht, hadden computers niet genoeg rekenkracht om deze uit te rekenen, maar inmiddels is dat wel het geval. Er is dus ook uitgebreid onderzoek gedaan om deze te vergelijken met de LIML, wat heel overzichtelijk is verzameld door Puhani (2000). Nawata (1994) vergelijkt de twee schatters met een Monte-Carlosimulatie die op dezelfde manier is opgezet als zijn onder-zoek van een jaar eerder. Hij komt tot de conclusie dat de FIML inderdaad beter is. Bij hoge ρ1, de situatie waar de LIML niet goed functioneert, werkt de FIML nog behoorlijk goed. Dat

blijkt vooral uit de variantie van de FIML-schatter, die veel kleiner is dan die van Heckmans tweestapssschatter.

De meeste andere wetenschappers vinden resultaten die bij de conclusie van Nawata aanslui-ten. Paarsch (1984) bijvoorbeeld concludeert ook dat de FIML in alle situaties beter is dan de Heckman-methode. Hij onderzoekt echter geen situaties waar er correlatie tussen de variabelen is, maar hij onderzoekt wat er gebeurt als de storingen niet normaal verdeeld zijn. Opvallend is

(8)

dat hij vindt dat ook als de storingen wel normaal verdeeld zijn, de tweestapsmethode slechter is. Dat komt doordat de Limited Information Maximum Likelihood slechts beperkte informatie gebruikt en daardoor niet efficiënt is. Leung en Yu (1996) komen echter tot een heel andere conclusie. Zij vinden geen duidelijke dominantie van de FIML. Daar is een goede verklaring voor: ze onderzoeken alleen situaties waar de ρ1gelijk is aan 0 of 0.5 en uit Nawata (1994) was

gebleken dat de dominantie van de FIML vooral tevoorschijn komt als de ρ1groter is dan 0.8.

Een heel andere schattingsmethode is het zogeheten Two Part Model (TPM). Het TPM is een alternatieve manier van kijken naar het sample-selection-probleem, die onder andere gebruikt en fel verdedigd wordt door Duan et al. (1983, 1984). Deze schatter komt in de rest van mijn onderzoek niet voor, maar hij wordt hier genoemd om te laten zien dat er ook alternatieven zijn die structureel anders naar het probleem kijken dan Heckmans tweestapmethode. Dit in tegen-stelling tot de FIML en de methoden van Newey en Coslett die hieronder besproken worden, die wel van hetzelfde onderliggende model uitgaan. Het model dat bij TPM geschat wordt, wordt als volgt omschreven. Hierbij is de notatie aangepast om zo veel mogelijk overeen te komen met die van de Heckman-methode.

y∗_1i|y∗_2i> 0 = x0_1iβ1+ u1i (14)

y∗_2i = x0_2iβ2+ u2i (15)

y_1i = y∗_1i als y∗_2i > 0 (16)

y_1i = 0 als y∗_2i ≤ 0 (17)

Het verschil van deze vergelijkingen met vergelijking (1) tot en met (4) is klein, maar de gevol-gen daarvan groot. De afhankelijke variabele wordt hier namelijk conditioneel gevol-genomen op een positieve waarde van y∗_2i. Dat zorgt ervoor dat de verwachte waarde wordt:

E(y∗_1i) = Φ(x0_2iβ2/σ2) ∗ (x1i0 β1) (18)

Dit is eenvoudiger te berekenen dan de Heckman-methode en het is minder afhankelijk van de normale verdeling. De Φ(..) kan namelijk eenvoudig vervangen worden door een andere verdeling als er het vermoeden bestaat dat de storingen niet normaal verdeeld zijn. Een nadeel is dat de interpretatie van de coëfficiënten moeilijker is. Het is een open vraag of men in de praktijk geïntereseerd is in de conditionele verwachting. Duan et al. (1984) vinden uiteraard

(9)

van wel, maar bijvoorbeeld Puhani (2000) is daar wat sceptischer over.

Deze twee alternatieven lijken dus op basis van de bovengenoemde onderzoeken beter te zijn dan de tweestapsmethode van Heckman. Ze hebben echter ook hun nadelen. De FIML blijft ingewikkelder om te berekenen en het verschil met de LIML is klein, aangezien hij voor een groot deel op dezelfde aannames berust. Het TPM heeft weer andere problemen, vooral op het gebied van interpretatie. Afgezien van deze twee schatters, zijn er nog meer alternatieven ontwikkeld. Ik noem hier alleen nog de methode van Cosslett (1991) en focus daarna op Newey (2009), aangezien dat de methode is die in de rest van dit onderzoek centraal staat.

2.4 De methode van Newey

De schatter van Newey (2009) is onder andere gebaseerd op die van Cosslett (1991). Beide methodes zijn een aanpassing van de Heckman-methode waarbij de Mills-ratio (zie vergelij-king (8)) vervangen wordt door een schatter die flexibeler is. Cosslett gebruikt een reeks van interval-dummy’s, waarbij op elk interval een aparte coëfficiënt geschat wordt. Deze methode valt echter niet binnen dit onderzoek.

Wat Newey doet, is de λ -functie schatten met behulp van een polynoom. Dat zorgt ervoor dat het uiteindelijk te schatten model er als volgt uitziet:

y_1i = x0_1iβ1+ K

∑

k=0 ηk∗ pk(x02iβ2) + ηi∗ (19) met η_i∗ = ∞

∑

k=K+1 η_k∗ p_k(x0_2iβ2) + ηi (20)

Hierbij zijn pk(..) bekende continue functies. Die functies kunnen in principe alles zijn, maar

er wordt ervoor gekozen om te nemen p_k(x) = [τ(x)]k−1, zodat het een polynoom is van orde k-1. Zo kan het effect op E(u1|u2> −x2i0 β2) consistent geschat worden, aangezien elke functie

benaderd kan worden door middel van een polynoom. Daarvoor is wel nodig dat K → ∞ als N→ ∞. Dat is ook de reden waarom de storingen uitgedrukt zijn als in vergelijking (20).

Voordat vergelijking (19) kan worden geschat, moet in de eerste stap β2 worden geschat.

Ook voor het schatten van die eerste stap stelt Newey (2009) een verbetering voor. Hij gebruikt namelijk een schatter die onafhankelijk is van de verdeling van de storingen uit het selectiecri-terium, zoals de methode van Klein en Spady (1993). Het valt echter buiten de opzet van dit onderzoek om die methode uit te leggen.

(10)

Het is echter nog niet gezegd wat voor functie τ(x) zou moeten zijn. Het is mogelijk om daar een eenvoudige lineaire functie te gebruiken, wat de makkelijkst te interpreteren resultaten oplevert. Newey (2009, p. 220) stelt verschillende functies voor, onder andere τ(x) = 2Φ(x) − 1 waarbij Φ(x) wederom de standaardnormale verdelingsfunctie is. Deze keuze wordt nagevolgd door onder andere Hussinger (2008, p. 735). Hier wordt dus niet aangenomen dat de storingen uit een normale verdeling afkomstig zijn, maar de geschatte resultaten x0_2iβˆ2 worden

gestrans-formeerd om het proces robuuster te maken.

Het grote voordeel van deze methode is dat het onderliggende verband benaderd kan worden, ook als dat niet lineair is dankzij het gebruik van de polynoom ∑∞

k=K+1ηk∗ pk(x02iβ2). Als K>1

is de vergelijking al niet lineair en zijn de grootste problemen met multicollineariteit dus al ontweken. Hier wordt ook al direct een nadeel zichtbaar. Hoe moet er bepaald worden hoe groot K is? Een goede optie daarvoor is om K steeds groter te kiezen, zolang de goodness of fit niet goed genoeg is. Newey (2009, p. S223) geeft verder geen suggesties welk criterium gebruikt kan worden om te bepalen wanneer K groot genoeg is.

Er is al eerder onderzoek gedaan naar deze methode. Hussinger (2008) doet dat op een tamelijk originele wijze. Zij onderzoekt, zoals in de inleiding vermeld, het effect van subsidies op innovatie. Daarvoor gebruikt ze vier verschillende methodes om zo tegelijk onderzoek te kunnen doen naar de verschillende resultaten van deze methoden. Ze maakt daarvoor gebruik van Heckman (1976, 1979) en van Newey (2009). Daarnaast gebruikt ze twee methodes die buiten dit onderzoek vallen, namelijk die van Cosslett (1991) en van Robinson (1984).

De resultaten die zij uit de verschillende methoden krijgt zijn helaas niet zo heel verschillend van elkaar. Alle vier de gebruikte methoden geven vergelijkbare resultaten en zelfs OLS wijkt nauwelijks af. Dat komt waarschijnlijk doordat er in de data minder selection bias zit dan van tevoren de verwachting was. Het is jammer dat door deze eigenschap van de data het vergelij-ken minder goed lukt. Daarom is het nuttig om een onderzoek te doen met gesimuleerde data, zodat de verschillen hopelijk wel zichtbaar zijn.

(11)

3 Onderzoeksopzet

Het doel van dit onderzoek is om te testen of de Newey-methode beter is dan de Heckman-methode. Dat gebeurt door middel van een computersimulatie. Voor deze Monte-Carlosimulatie wordt van het volgende datagenererende proces gebruik gemaakt:

y∗_1i = β0+ β1x1i+ u1i (21)

y∗_2i = α0+ α1x2i+ u2i (22)

y_1i = y∗_1i als y∗_2i> 0 (23)

y_1i = 0 als y∗_2i≤ 0 (24)

Daarbij zijn de achterliggende verdelingen als volgt:   x1 x₂   ∼ BivUni f     10 10  ,   100/3 ρ1∗ 100/3 ρ1/ ∗ 100/3 100/3     (25)   u1 u₂   ∼ BN     0 0  ,   25 ρ2∗ 25 ρ2∗ 25 25     (26)

Dat betekent dat de storingen normaal verdeeld zijn met verwachting 0 en standaardafwijking 5 en correlatie tussen u1en u2gelijk aan ρ2. Voor x1en x2betekent het dat ze uniform verdeeld

zijn tussen 0 en 20 en dat hun correlatie gelijk is aan ρ1.

Hierbij wordt α0zodanig gekozen dat die gelijk is aan −E(x2iβ2) = −10, zodat de verwachte

waarde van y2i gelijk is aan 0, zodat de verwachte truncatie op 0.5 wordt gesteld. Er is al

onderzoek gedaan naar wat er gebeurt als de truncatie verandert (bijvoorbeeld door Paarsch, 1984), maar dat leverde geen interessante resultaten op. Om die reden wordt dat hier niet verder onderzocht.

De echte waarde van β0wordt gesteld op 10. Bij Nawata (1993, 1994) wordt -10 gebruikt,

maar dat zorgt ervoor dat er ook veel negatieve waarden van y1 voorkomen. Dat is niet in

overeenstemming met de situaties uit de praktijk waarvoor de methodes gebruikt worden, zoals lonen of subsidies. De echte waarden van β1en β2zijn zonder verdere achterliggende gedachte

op 2 en 1 gesteld.

Voor het schatten van de eerste stap wordt bij beide methodes Probit gebruikt. Aangezien al-leen storingen afkomstig uit een normale verdeling worden gebruikt, is Probit de beste methode. Als er echter ook storingen uit een andere verdeling onderzocht zouden worden, zou er gebruik

(12)

moeten worden gemaakt van semiparametrische methoden. Dat gebeurt in deze simulatie echter niet.

Wat hier wel wordt onderzocht is hoe goed de methoden van Newey en Heckman toepasbaar zijn onder verschillende omstandigheden. De situaties die onderzocht worden zijn:

1. Verschillende waarden voor ρ1, namelijk 0, 0.25, 0.5, 0.8, 0.9, 0.95

2. Verschillende waarden voor ρ2, namelijk 0, 0.5, 0.9

De focus ligt vooral op punt 1, aangezien Nawata (1993) al heeft aangetoond dat de schat-tingen met de Heckman-methode bij grote waarden van ρ1 slecht zijn. De verwachting is dat

Newey op dit vlak beter presteert.

Hoe goed de beide methodes werken, wordt bepaald op basis van de grootte van de bias van de schatters en op basis van de variantie van de schatters. De datagenererende waarden van de steekproefgrootte en het aantal replicaties, ten slotte, zijn als volgt:

replicaties = 500, steek proe f grootte = 200 (27)

De omvang van de steekproef is relatief klein, wat betekent dat de graad van de polynoom bij Newey, K, ook klein moet zijn. In dit onderzoek is ervoor gekozen om die op 3 te zetten. De verwachting is namelijk dat bij grotere steekproeven en dus grotere K de polynoom sowieso een goede benadering geeft. Dat is bij een kleine waarde zoals deze minder waarschijnlijk en dus interessanter om te onderzoeken.

4 Resutlaten

Op de volgende pagina’s staan de tabellen met daarin de resultaten van de computersimulatie. Het zijn drie tabellen, één voor elke waarde van ρ2 die is onderzocht. Voor elke waarde van

ρ1 zijn in alle drie de tabellen de schattingen van de coëfficiënten voor β0 en β1 met beide

schattingsmethoden opgenomen. Bij elk van die schattingen staat vervolgens de waarde van de schatter, gemiddeld over de vijfhonderd replicaties. Vervolgens zijn ook de standaardafwijking, mediaan en het 25ste en 75ste kwantiel daarvan weergegeven. De standaardafwijking in de tabellen is dus niet de geschatte waarde van de populatiestandaardafwijking, maar de afwijkijng van de schatting in de verschillende replicaties.

Het eerste wat opvalt in tabel 1, is dat alle schattingen zuiver lijken. De gemiddelde schattin-gen van β1wijken bij alle waarden van ρ1niet meer dan 1 procent af van de echte waarde die 2

(13)

Tabel 1: Waarden voor ρ2= 0

ρ1 Coëfficiënt Gemiddelde Standaardafwijking Mediaan 25ste Kwantiel 75ste Kwantiel

0.0 Heckman β0 10.04928 1.142146 10.01048 9.222997 10.81447 β1 1.994641 0.09206999 1.993431 1.937334 2.05622 Newey β0 10.05016 1.286621 9.991709 9.133455 10.88827 β1 1.995134 0.09386301 1.993641 1.930415 2.05832 0.25 Heckman β0 9.984764 1.372296 9.990111 9.089914 10.82684 β1 1.999684 0.09508866 2.003 1.936822 2.066302 Newey β0 9.966453 1.340011 9.929256 9.047338 10.94342 β1 2.000457 0.09524712 2.003192 1.938396 2.064287 0.5 Heckman β0 9.982441 1.639929 10.00024 8.822426 11.07079 β1 1.99954 0.10514855 1.995379 1.93343 2.070157 Newey β0 9.955659 1.429426 9.91564 9.033801 10.98346 β1 2.000139 0.1070892 1.994636 1.934542 2.072339 0.8 Heckman β0 10.03243 2.439909 10.06722 8.370699 11.64873 β1 1.99671 0.14338298 1.991109 1.907535 2.088258 Newey β0 9.99287 1.797852 9.905858 8.808512 11.14733 β1 1.997452 0.1520709 1.999319 1.904889 2.087829 0.9 Heckman β0 10.02919 3.35987 10.09388 7.719075 12.13427 β1 1.997406 0.18998589 1.983853 1.862476 2.136377 Newey β0 10.00347 2.31596 9.824477 8.54265 11.52369 β1 1.99826 0.2078577 1.993846 1.86334 2.140266 0.95 Heckman β0 10.0675 4.443537 10.22915 6.984725 12.94142 β1 1.995358 0.24543771 1.992248 1.829995 2.172602 Newey β0 10.01275 3.016777 9.963536 8.041399 11.97708 β1 1.997308 0.2838301 1.9959 1.80841 2.191568

naarmate ρ1 groter wordt. Dat blijkt uit het feit dat de standaardafwijking steeds groter wordt

en tevens daaruit dat de kwartielen steeds verder uit elkaar komen te liggen.

Wat verder opvalt, is dat de verschillen tussen de Heckman-methode en de Newey-methode uiterst klein zijn. Het lijkt erop dat de standaardafwijking bij Newey iets groter is, maar dat kan ook komen doordat bij het berekenen van deze waarden drie termen zijn opgenomen voor de

(14)

Tabel 2: Waarden voor ρ2= 0.5

polynoom van Newey. Dat lijkt weliswaar weinig, maar in de meeste gevallen was alleen de eerste coëfficiënt van de polynoom significant, dus K was wellicht zelfs te hoog voor de grootte van de steekproef. Dat zou betekenen dat er sprake is van enige mate van overfitting, wat ervoor zorgt dat de schattingen onnauwkeuriger worden en daardoor de standaardafwijkingen hoger.

(15)

Tabel 3: Waarden voor ρ2= 0.9

Het is echter vreemd dat bij de schattingen van β0 het tegenovergestelde aan de hand lijkt te

zijn. Ook deze waarden lijken zuiver voor beide schattingsmethoden en ook bij deze schattin-gen neemt de variantie toe als de correlatie tussen x1en x2toeneemt. Het grote verschil is echter

dat hier de standaardafwijking bij de Heckman-methode groter is dan bij de Newey-methode. Dat impliceert dat Newey betrouwbaarder resultaten voor de constante schat. Bij Heckman is

(16)

namelijk als ρ1= 0.9 in bijna 25% van de gevallen de schatter hoger dan 13 en in 25% van de

gevallen kleiner dan 7. Bij Newey is dat 12 en 8, wat aanzienlijk beter is.

Als ρ2 gelijk is aan 0.5 of 0.9 ligt de situatie heel anders. De schattingen van β1 zijn heel

vergelijkbaar met de situatie van ρ2= 0, waarbij nog opgemerkt kan worden dat de

standaard-afwijking iets kleiner is. Dat is ook logisch, want als de storingen gecorreleerd zijn, is er in totaal dus iets minder storing. Dan kan er dus betrouwbaarder geschat worden. De schatting van β0 echter, wijkt bij Newey af van de populatiewaarde van 10. Zeker als ρ2= 0.9 wijst

het gemiddelde van 13.5 op een systematische bias. Dat komt doordat er een polynoom wordt toegevoegd aan de regressie in de tweede stap. Die polynoom heeft zelf ook een constante, wat als gevolg heeft dat het onmogelijk wordt om altijd consistente schattingen voor β0te krijgen.

Dit betekent dat de geschatte waarde van β0 zijn interpretatie verliest. De schattingen bij

Heckman zouden deze onzuiverheid niet moeten hebben en dat blijkt ook uit de tabellen. Het lijkt erop dat bij hoge waarden van ρ1de constante gemiddeld iets te laag geschat wordt. Een

afwijking van 0.12 bij een standaardafwijking van 3.48, zoals bij ρ2 = 0.9 en ρ1= 0.95, is

echter niet significant. Dat is verwonderlijk, want uit eerdere onderzoeken (zoals Nawata, 1993, 1994) was gebleken dat de schattingsresultaten met Heckmans tweestapsmethode significant verschillen van de echte waarde, als de correlatie tussen de verklarende variabelen groot is. Dat dat hier niet het geval is, wordt mogelijk veroorzaakt doordat de variantie van de storingen niet bijzonder groot is, waardoor de afwijking te klein is om waar te nemen.

Het belangrijkste resultaat dat uit de tabellen naar voren komt, is dat de verschillen tussen de methoden in de onderzochte situaties bijzonder klein zijn. Het is blijkbaar mogelijk om zelfs met een polynoom van een lage graad schattingen te krijgen voor β1die even goed zijn als met

de inverse Mills-ratio. In alle onderzochte situaties waren de storingen echter normaal verdeeld. Dat is voor beide methoden gunstig, aangezien dan in de eerste stap Probit gebruikt kan wor-den. Voor de tweestapsmethode van Heckman is het echter een noodzakelijke voorwaarde om consistente resultaten te krijgen, terwijl dat bij Newey niet het geval is. De methoden zijn dus even goed toepasbaar in deze situatie, waar Newey in het nadeel is door de normale verdeling en door de lage graad van de polynoom. Hierdoor is de verwachting dat in andere situaties Newey wel betere resultaten oplevert. Verder onderzoek zal dit moeten uitwijzen.

(17)

5 Conclusie

Het doel van dit onderzoek was om te kijken welke methode de beste oplossing geeft voor sample selection bias. Die bias komt voort uit het feit dat in bepaalde economische situaties een deel van de data niet kan worden waargenomen, waardoor technieken die daar niet expliciet rekening mee houden een vertekend beeld geven. De gebruikelijke oplossing is om daarvoor de Heckman tweestapsmethode te gebruiken. Er is echter een aantal theoretische bezwaren tegen de methode van Heckman. Deze methode is namelijk alleen bruikbaar met normaal verdeelde storingstermen en de verwachting is dat zij problemen ondervindt met multicollineariteit, als er grote correlatie is tussen de verklarende variabelen uit de selectievergelijking en de verklarende variabelen. Ook uit eerdere Monte-Carlosimulaties is gebleken dat de methode van Heckman niet altijd optimaal presteert.

Daarom zijn er verschillende alternatieve methoden ontwikkeld, waaronder die van Newey. Deze methode maakt gebruik van een polynoom om de onzuiverheid die voortkomt uit de cen-surering van de data, te schatten. In dit onderzoek werden door middel van computersimulaties de verschillen tussen de resultaten die de methode van Heckman en de methode van Newey opleveren, met elkaar vergeleken. Daarbij lag de focus op wat er zou gebeuren als de correlatie tussen de storingen en tussen de verklarende variabelen groot zijn.

Verrassend genoeg is gebleken dat de Heckman tweestapsmethode niet onzuiver werd, ter-wijl dat in eerdere onderzoeken wel gebeurde. Verder bleek de methode van Newey niet in staat om de constante te schatten als de correlatie tusssen de storingen hoog is. Als een onderzoeker dus geïnteresseerd is in de echte waarde van de constante, is het beter om de Heckman-methode te gebruiken dan de Newey-methode. Daar staat tegenover dat de twee methoden met betrek-king tot β1resultaten geven die nauwelijks van elkaar verschillen.

Het lijkt er dus op dat de methode van Newey niet beter is dan de tweestapsschatter van Heckman. Lang niet alle situaties zijn echter hier onderzocht. In alle onderzochte gevallen waren de storingen afkomstig uit een normale verdeling. Bij andere verdelingen zullen de re-sultaten zeer waarschijnlijk anders zijn. Ook is het aantal termen dat in de Newey-polynoom is gebruikt in alle situaties hetzelfde wat ook invloed heeft op de resultaten. De verwachting is dat in andere gevallen de verschillen tussen de methoden beter zichtbaar zijn. Vooralsnog echter geven beide methodes even goede resultaten.

(18)

6 Bibliografie

Coslett, S.R. (1991). Semiparametric estimation of a regression model with sample selecti-vity. In Nonparametric and semiparametric estimation methods in econometrics and statistics, Barnett, W.A., Powell, J. & Tauchen, G.(eds.), Cambridge(U.K.): Cambridge University Press, 175-197.

Duan, N., Manning, W.G., Morris, C.N. & Newhouse, J.P. (1983). A Comparison of Alter-native Models for the Demand for Medical Care. Journal of Business & Economic Statistics, 1(2), 115-126.

Duan, N., Manning, W.G., Morris, C.N. & Newhouse, J.P. (1984). Choosing Between the Sample-Selection Model and the Multi-Part Model. Journal of Business & Economic Statistics, 2(3), 283-289.

Heckman, J.J. (1976). The Common Structure of Statistical Models of Truncation, Sample Selection and Limited Dependent Variables and a Simple Estimator for Such Models. Annals of Economic Social Measurement, 5(4), 475-492.

Heckman, J.J. (1979). Sample Selection Bias as a Specification Error. Econometrica, 47(1), 53-161.

Heij, C., de Boer, P., Franses, P.H., Kloek, T. & van Dijk, H.K. (2004). Econometric Methods with Applications in Business and Economics, New York: Oxford University Press Inc.

Hussinger, K. (2008) R&D and subsidies at the form level: An application of parametric and semipara-metric two-step selection models. Journal of Applied Econometrics 23(1), 729-747.

Klein, R. W., Spady, R. H. (1993) An Efficient Semiparametric Estimator for Binary Res-ponse Models. Econometrica 61 (2), 387-421.

Leung, S.F., Yu S. (1996). On the Choice Between Sample Selection and Two-Part Models. Journal of Econometrics, 72(1), 197-229.

Little, R.J.A. & Rubin, D.B. (1987). Statistical Analysis with Missing Data. New York: John Wiley & Sons.

Nawata, K. (1993). A Note on the Estimation of Models with Sample Selection Biases. Economics Letters, 42(1), 15-24.

Nawata, K. (1994). Estimation of Sample Selection Bias Models by the Maximum Likeli-hood Estimator and Heckman’s Two-Step Estimator. Economics Letters, 45(1), 33-40.

(19)

Econo-metrics Journal, 12(1), S217-S229.

Paarsch, H.J. (1984). A Monte Carlo Comparison of Estimators for Censored Regression Models. Journal of Econometrics, 24(1), 197-213.

Puhani, P.A. (2000). The Heckman Sorrection for Sample Selection and its Critique. Journal of Economic surveys, 14(1), 53-68.

Robinson, P. (1988). Root-N-consistent sempiparametric regression. Econometrica 56(4), 931-954.