-
Notifications
You must be signed in to change notification settings - Fork 0
/
weakiv.ado
7703 lines (6993 loc) · 261 KB
/
weakiv.ado
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
*! weakiv 2.4.07 1Oct2015
*! authors Finlay-Magnusson-Schaffer
* Usage:
* <eqn>, <options> = estimate equation as model=linear/ivprobit/ivtobit.
* <nothing>, <options> = uses model estimated by ivregress/ivreg2/etc. in memory
* if no such model in memory, works as Stata replay and
* replays last weakiv estimates. useful for graph tweaking.
* <nothing>, estusewald(*) = same as above but uses stored model
* <nothing>, version = reports version number, does nothing else
* Version notes:
* 1.0.01 31Jul2013. First complete working version
* 1.0.02 04Aug2013. Bug fix in kwt (was ignoring option)
* 1.0.03 05Aug2013. Renamed weakiv with modified syntax. Minor program restructuring.
* 1.0.04 07Aug2013. J cset reported. Open-ended and entire-grid csets noted in output.
* Hyperlinks in output table added; link to sections of help file.
* "all" option added for graph(.).
* 1.0.05 11Aug2013. Added support for ivreg2h (Lewbel-type generated IVs).
* Fixed bug in combination of replay() and graph(all).
* 1.0.06 22Sep2013. Major update. Added support for 2-endog-regressor case.
* Graphing for K=2 case requires Stata 12 (contour) and -surface-.
* Switch from use of ranktest to use of avar.
* level option for graphing now allows multiple levels;
* with replay can change level vs. estimation table.
* Various coding reorganizations, fixes and tidying up.
* Bug fix for iid closed-form CIs when only included exog is constant.
* Minor bug fix for CLR - numerical precision issues mean it can be
* very small and <0, so now take abs(.)
* 1.0.07 11Nov2013. Minor bug fix for K=2 (would crash in grid search if exactly IDed)
* 1.0.08 15Jan2014. Added support for FE and FD estimation by xtivreg2 and xtivreg.
* Included recoding for temporary variables to support TS operators.
* Added number of observations to table output.
* Added contourline (default) and contourshade to contour options
* Added kjlevel(.), arlevel(.), jlevel(.) options.
* Saved level macros now level and levellist (if >1 specified)
* Behavior of level(.) option is now that first (not highest) is used for tests.
* Added check for Stata-legal confidence levels (>=10, <=99.99)
* Save cluster numbers and variables as e(.) macros
* Add support for Kleibergen LM method, both iid and non-iid cases. Implies partialling-out.
* NB: Kleibergen non-iid formulae do not reduce to iid formula when iid covariances used.
* Fixed bug in closed-form CIs with dofminus
* 1.0.09 16Feb2014 Fixed bug in user-set levels for J and AR tests.
* Fixed display of excessive precision for kwt.
* 1.0.10 23Mar2014 Added strong(.) option. Various code tweaks.
* Fixed minor contouronly/surfaceonly bug (ignored in interactive mode)
* 1.0.11 23May2014 Fixed contourshade bugs; can now control starting/ending color.
* Fixed minor contourlevel bugs (was ignoring additional specified levels)
* 2.0.00 26Jun2014 Major internal rewrite to accommodate xtabond2 and aribtrary number of endogenous vars.
* Changes include: always preserve/restore; always partial-out in linear models;
* support for factor variables; use of sclass programs for parsing;
* promoted to requiring Stata 11.2 (because of _rmcollright bug under version control v=10.1);
* tests reported for arbitrary number of weakly and strongly identified coeffs;
* CIs and graphs available for #weak=1 and #weak=2 and arbitrary number of strongly IDed;
* specification of null options changed to single numlist;
* misc minor changes to saved results
* 2.0.01 28Jun2014 Fix in xtabond2 parsing code to catch whether data in levels, diffs or both are used.
* xtabond2 parsing catches whether there are zero endogenous regressors.
* Fixed bug in npd flag not being set to 1 if npd matrix encountered in test of specific null.
* Fixed bug in use of marksample.
* 2.0.02 30Jun2014 Bug fixes: levels vs. diffs data with xtabond2; whether xtabond2 has K>0; output formatting;
* saved matrices to enable xtabond2 to work as postestimation command.
* Report various ob and var counts with xtabond2. Report list of strongly identified.
* Added EQxtabond2(.) option to override default of diff/level/sys data.
* Added md option. Default method now LM for linear models, MD for others.
* 2.0.03 4July2014 Errors in pre-weakiv estimation command now capture-d rather than quietly-ed.
* Stata's built-in _rmcollright provided with Stata 12.1 (version 1.3.0 07apr2009) has
* a bug that doesn't allow it to accept weights, so replaced it with weakiv_rmcollright
* that has same functionality but with bug fixed.
* Mikusheva-Poi CLR code revised to assuming partialling-out.
* Fixed bug that arose if the weighting variable was also a regressor.
* Added support for xtabond2 with weights and with classical VCE.
* Added support for xtabond2 with non-panel-id clustering variable (version forthcoming)
* 2.0.04 5July2014 Fixed bug in reported names of endog/inexog regressors with xtabond2 and eq(.) option meant
* subset of data (just diff or lev out of total sys estimation) used.
* 2.0.05 5July2014 Caught minor bug where weakiv would crash with xtabond2 if #endog=0
* 2.0.06 17July2014 Catch and warn if graphs requested when #endog>2. Added g_min, g_max, g_avg macros.
* Recoded CI table in Mata to avoid Stata max matsize limit. Total #points in output.
* 2.1.00 28July2014 Rewrite of code to allow grid searches in arbitrary number of dimensions.
* Rewrite of code to report projection-based inference
* Changes of options: gridpoints, gridlist, gridmin, gridmax
* Macro name change to grid_desc; new macro added points_desc
* Fixed bug in column order of stored nulls
* Fixed bug in reporting of CIs if interval started on last gridpoint
* 2.2.00 5Aug2014 Added subset AR test. Added projection-based inference. More code rewrites.
* Added cuestrong option for LIML and CUE estimation for strongly-IDed coeffs.
* Fixed bug in CIs with LM option - was using MD. Strongly-IDed coeffs now saved in citable.
* 2.2.01 7Aug2014 Added scatter option for 2-way CI sets. Added 2-way projection inference option.
* Projection-based CIs saved as e(.) matrices.
* 2.3.00 11Aug2014 rk statistic (K-P rank test) now always available. CLR stat now always available.
* Added p-values for rk statistic.
* p-values for CLR stat available through simulation.
* Fixed bug in replay with project(.)
* Rationalized code (single routine for computing tests).
* Contourshade now default for 2-D graphs.
* 2.3.01 21Aug2014 Added testexog(.) option for including exogenous regressors in tests
* 2.3.02 24Aug2014 Fixed bug: 2-way graphs were leaving behind temporary variable
* Added testid option to report tests of underidentification
* 2.3.03 9Sept2014 Added cuepoint option for reporting CUE point estimates
* Added CUE-MD method for CUE with Wald/MD tests
* 2.3.04 12Sept2014 Fixed bug where project2(.) didn't trigger grid construction
* cuepoint option triggers centering of grid on CUE rather than original Wald estimator
* Allow spec with one gridpoint in a dimension
* (NB: first public release of version 2)
* 2.3.05 2Oct2014 Fixed bug in reporting of #exexog with xtabond2 - was reporting 0
* Added saved macro for #exexog
* 2.4.00 1Jan2015 Fixed bug in unsaved macro with list of CI tests.
* Changed replay to use saved macro list of tests.
* Added Wald to list of projection-based tests.
* Minor recoding to separate Wald CI code from CI code for rest of tests.
* Rewrite of code for extracting CIs; rejection indicators no longer saved;
* replay mode now triggers recalculation of CIs and rejection indicators.
* Minor changes to names of saved matrices (betas and VCEs)
* Rewrite of code for storing and working with CI tables; now only exist in Mata
* until end of program, when stored as e(.) if #rows <= 32,767.
* If citable is too large to start as e(citable),
* split into separate e(.) matrices of 10k rows each.
* 2.4.01 18Jan2015 Changed call to (xt)ivreg2 to use nooutput option so that warning messages are displayed.
* Checks whether ivreg2h is called with panel-data option (not currently supported)
* 2.4.02 9Feb2015 Fixed bug in gridlist option; was ignoring supplied list.
* Changed cuepoint option so that it triggers inclusion of CUE beta in grid
* 2.4.03 11Feb2015 Added support for updated xtabond2 - e(clustid) became e(clustid1)
* Fixed bugs causing tokenize to choke if string had commas, e.g., "l(1,2).somevar"
* 2.4.04 26Mar2015 Added support for ivreg2 with factor variables; also general recoding of support
* for factor vars using -fvstrip- utility.
* Fixed bug in xtabond2 parsing code; now allows possibly-empty inexog list
* Fixed small bug in xtabond2 - wasn't reporting correct number of clusters in rare cases
* Fixed reporting of stats when encountering NPDs - negative stats are changed to missing
* 2.4.05 15Jun2015 Update to internal utility fvstrip
* 2.4.06 1July2015 Bug fix - projection intervals were not accepting non-integer K-J test level.
* Code tweak - "capture drop __wt" changed to "capture drop __wt1"
* 2.4.07 1Oct2015 Rewrite of code to parse xtabond2.
* 2.4.08 16Feb2016 Bug fix - cuepoint option with exactly-IDed model was causing LIML code to crash.
* to do: note possible upper limit to points that contour can render; use scatter option instead
* use mat colnames from saved matrices with xtabond2
* recommend using LIML or CUE with cuepoint so that automatic graph grid is centred correctly.
* fix bug that disallows using a string variable for the cluster option
* consider replacing cholesky(invsym(psi)) with rspi_inv = matpowersym(psi, -0.5)
program define weakiv, eclass byable(recall) sortpreserve
version 11.2
local lversion 02.4.08
local avarversion 01.0.04
local ranktestversion 01.3.03
checkversion_avar `avarversion' // Confirm avar is installed (necessary component).
if replay() { // No model provided before ",", but possibly options
syntax [, ///
VERsion ESTUSEwald(name) ///
* ///
]
if "`version'" != "" { // Report program version number, then exit.
di in gr "`lversion'"
ereturn clear
ereturn local version `lversion'
exit
}
if "`e(cmd)'"=="weakiv" & /// If previous e(cmd) was weakiv, and Wald model is NOT provided,
"`estusewald'" == "" { // then replay last weakiv results with any new options, then exit.
weakiv_replay `0' // Useful for tweaking graph options.
exit // weakiv_replay syntax will catch illegals.
}
if "`estusewald'" != "" { // Wald model is provided by user in estusewald(.) option,
est restore `estusewald' // so make current and then proceed to weak-IV-robust estimation.
}
// If Wald model is NOT provided then it must already be in memory
// as e(.) results, so proceed to weak-IV-robust estimation.
}
else { // Estimate model specified by user before ",".
estimate_model `0' // Estimate model, then proceed to weak-IV-robust estimation.
} // end replay()
************************* ASSEMBLE MODEL AND OPTION SPECS *****************************************
* Model in memory is now main model for Wald tests etc.
* Weak IV-robust estimation inherits same characteristics (robust etc.).
* Save command line
local cmdline "weakiv `0'"
* esample = touse, except for xtabond2 where esample = original e(sample) in estimation
* use touse unless need to refer explicitly to original sample (e.g. when expanding factor vars)
tempvar esample
qui gen byte `esample'=e(sample)
* If data already preserved, OK; if preserve fails for other reasons, exit with error
capture preserve
if _rc > 0 & _rc~=621 {
di as err "Internal weakiv error - preserve failed"
exit 498
}
* Create temp variables that will be changed by subroutines
tempvar touse wvar clustvar1_t clustvar2_t
qui gen byte `touse'=.
qui gen double `wvar'=.
qui gen long `clustvar1_t'=.
qui gen long `clustvar2_t'=.
* Get weakiv options from command line
* Pass additional args to get_option_specs; bind needed in case of constructs like "l(1,2).abmi"
gettoken first opts : 0 , parse(",") bind
if "`first'"~="," { // args in macro `first' so run again
gettoken first opts : opts , parse(",") bind // to get rid of comma
}
* Get model specs from active model
* Clear any extraneous sreturn macros first
sreturn clear
get_model_specs, /// Gets model specs from prev model. Catches some illegals.
`opts' /// Will need to pick up strong(.) from options and check vs. endo
touse(`touse') /// Applies to possibly-expanded current estimation sample
esample(`esample') /// Applies to original data estimation sample
wvar(`wvar') ///
clustvar1_t(`clustvar1_t') ///
clustvar2_t(`clustvar2_t')
local waldcmd "`s(waldcmd)'" // command for original Wald (non-weak-ID-robust) estimation
local model "`s(model)'"
local xtmodel "`s(xtmodel)'" // empty if not estimated using panel data estimator
local ivtitle "`s(ivtitle)'"
local depvar "`s(depvar)'"
local depvar_t "`s(depvar_t)'"
local inexog "`s(inexog)'"
local inexog_t "`s(inexog_t)'"
local tinexog "`s(tinexog)'"
local tinexog_t "`s(tinexog_t)'"
local endo "`s(endo)'"
local endo_t "`s(endo_t)'"
local wendo "`s(wendo)'"
local wendo_t "`s(wendo_t)'"
local sendo "`s(sendo)'"
local sendo_t "`s(sendo_t)'"
local csendo "`s(csendo)'"
local csendo_t "`s(csendo_t)'"
local pwendo "`s(pwendo)'" // Projection-based inference. Don't need _t version.
local pwendo2 "`s(pwendo2)'"
local csendo "`s(csendo)'" // Not in subset AR test(complement to set). Don't need _t version.
local exexog "`s(exexog)'" // Empty if xtabond2
local exexog_t "`s(exexog_t)'"
local gmminsts1 "`s(gmminsts1)'" // Specific to xtabond2
local gmminsts2 "`s(gmminsts2)'" // Specific to xtabond2
local ivinsts1 "`s(ivinsts1)'" // Specific to xtabond2
local ivinsts2 "`s(ivinsts2)'" // Specific to xtabond2
local noconstant "`s(noconstant)'" // "noconstant" if no constant in original model specification
local cons "`s(cons)'" // 0 if no constant in original model
local nendog "`s(nendog)'" // These counts will be overwritten after collinears etc. are removed
local nsendog "`s(nsendog)'" // #strongly identified endog; also boolean for strongly-IDed
local nwendog "`s(nwendog)'" // #weakly identified endog
local npwendog "`s(npwendog)'" // #wendo to construct project-based CIs; also boolean for projection method
local ncsendog "`s(ncsendog)'" // #wendo in complement to subset test; also boolean for subset test
local nexexog "`s(nexexog)'"
local ninexog "`s(ninexog)'"
local ntinexog "`s(ntinexog)'"
local nexog "`s(nexog)'"
local small "`s(small)'"
local robust "`s(robust)'"
local cluster "`s(cluster)'" // = cluster( <varname> ) or cluster( <varlist> )
local clustvar "`s(clustvar)'" // ="clustvar1" if 1-way, ="clustvar1 clustvar2" if 2-way
local clustvar1 "`s(clustvar1)'" // clustvar1_t already defined
local clustvar2 "`s(clustvar2)'" // clustvar2_t already defined
local bw "`s(bw)'"
local kernel "`s(kernel)'"
local llopt "`s(llopt)'"
local ulopt "`s(ulopt)'"
local asis "`s(asis)'"
local wtexp "`s(wtexp)'"
local wtype "`s(wtype)'"
local exp "`s(exp)'"
local wf "`s(wf)'"
local N "`s(N)'"
local N_clust "`s(N_clust)'" // With 2-way clustering, is min(N_clust1,N_clust2)
local N_clust1 "`s(N_clust1)'"
local N_clust2 "`s(N_clust2)'"
local N_g "`s(N_g)'" // #panel groups
local g_min "`s(g_min)'"
local g_max "`s(g_max)'"
local g_avg "`s(g_avg)'"
local singleton "`s(singleton)'" // #panel singletons
local dofminus "`s(dofminus)'" // 0 unless set by ivreg2
local psd "`s(psd)'"
local vceopt "`s(vceopt)'"
local note1 "`s(note1)'"
local note2 "`s(note2)'"
local note3 "`s(note3)'"
local iid "`s(iid)'"
************************* Prep for transformations **************************************
* If TS or FV operators used, replace with temporary variables.
* Ensures that any transformation are applied to temporary vars.
* Must do at this level since if done in subroutines, temp vars would disappear after return
* Extension "_t" is list with temp vars.
* xtabond2 _t varlists are already transformed
* other _t lists are original varlists and will be overwritten
* Varlists come with full set of factor variables
* fvrevar will create a full set, so some can be all zeros,
* including the base variable.
* Call fvrevar 1 variable at a time so that repetitions lead to same tempvar.
* varlist will have full original set of variables
* varlist1 will have variables with collinear and zeroed-out vars dropped
* varlist_t will have varlist1 with tempvar equivalents
* nvarlist will have number of variables after collinears dropped
* Process depvar
fvrevar `depvar_t', substitute
local depvar_t "`r(varlist)'"
* Process endogenous
local templist
forvalues i=1/`nendog' {
local vn : word `i' of `endo_t'
fvrevar `vn', substitute
local templist "`templist' `r(varlist)'"
}
local endo_t : list clean templist // new varlist with tempvars
* Now remove collinears and zeroed-out factor variables etc.
clean_varlist `wtexp' if `touse', vlist(`endo') vlist_t(`endo_t') `noconstant'
local endo_c "`r(vlist_c)'" // cleaned original vlist with collinears removed etc.
local endo_t "`r(vlist_c_t)'" // corresponding temp vars
local nendog =r(nvars_c) // update count
* Process excluded exogenous
local templist
forvalues i=1/`nexexog' {
local vn : word `i' of `exexog_t'
fvrevar `vn', substitute
local templist "`templist' `r(varlist)'"
}
local exexog_t : list clean templist // new varlist with tempvars
* Now remove collinears and zeroed-out factor variables etc.
clean_varlist `wtexp' if `touse', vlist(`exexog') vlist_t(`exexog_t') `noconstant'
local exexog_t "`r(vlist_c_t)'"
local nexexog =r(nvars_c) // update count
* Process included exogenous
local templist
forvalues i=1/`ninexog' {
local vn : word `i' of `inexog_t'
fvrevar `vn', substitute
local templist "`templist' `r(varlist)'"
}
local inexog_t : list clean templist // new varlist exexog_t with tempvars
* Now remove collinears and zeroed-out factor variables etc.
clean_varlist `wtexp' if `touse', vlist(`inexog') vlist_t(`inexog_t') `noconstant'
local inexog_t "`r(vlist_c_t)'"
local ninexog =r(nvars_c)+`cons' // update count
* Now overid will be correct
local overid = `nexexog' - `nendog'
* Finally, varlists that appear via options. tempvars will already exist.
local templist
foreach var in `wendo_t' {
fvrevar `var', substitute
local templist "`templist' `r(varlist)'"
}
local wendo_t : list clean templist
local templist
foreach var in `sendo_t' {
fvrevar `var', substitute
local templist "`templist' `r(varlist)'"
}
local sendo_t : list clean templist
local templist
foreach var in `csendo_t' {
fvrevar `var', substitute
local templist "`templist' `r(varlist)'"
}
local csendo_t : list clean templist
local nwendog : word count `wendo_t' // update count
local nsendog : word count `sendo_t' // update count
local ncsendog : word count `csendo_t' // update count
local nexog = `nexexog'+`ninexog' // update count
* If dropped from endo, drop from wendo and sendo as well
local todrop : list endo - endo_c
local wendo : list wendo - todrop
local sendo : list sendo - todrop
local csendo : list csendo - todrop
***************************** WALD MODEL **********************************************
* Will use later for Wald tests and delete unless eststorewald specified.
* Also used for grid options
tempname waldmodel fullbeta var_fullbeta ebeta wbeta var_wbeta var1 var2
* For Wald tests and construction of stats and grids.
* Extract #nendog subvector of coeffs and corresp submatrix of VCV:
local cnames : colnames(e(b)) // not colfullnames since ivtobit has eqn names we don't want
local cnames : subinstr local cnames "bn." ".", all // strip out base notation
local cnames : subinstr local cnames "b." ".", all
local cnames : subinstr local cnames "o." ".", all
local cnames : list clean cnames // list of all regressors
foreach vn of local endo {
local pos : list posof `"`vn'"' in cnames
local endopos "`endopos' `pos'"
}
local endopos : list clean endopos // list of columns in which all endog coeffs appear
foreach vn of local wendo {
local pos : list posof `"`vn'"' in cnames
local wendopos "`wendopos' `pos'"
}
local wendopos : list clean wendopos // list of columns in which weakly-IDed coeffs appear
* Assemble parameter vector and VCE for weakly-identified coeffs
* ebeta is full set of Wald estimates; wbeta is weakly-ID only
mat `fullbeta' = e(b) // fullbeta has all coeffs including exogenous
mat `var_fullbeta'=e(V)
foreach cn of local endopos { // ebeta has all endog coeffs
mat `ebeta' = nullmat(`ebeta') , `fullbeta'[1,`cn']
}
mat colnames `ebeta' = `endo'
foreach cn of local wendopos { // wbeta has only weakly-ID coeffs
mat `wbeta' = nullmat(`wbeta') , `fullbeta'[1,`cn']
mat `var1' = nullmat(`var1') , `var_fullbeta'[1...,`cn']
}
foreach cn of local wendopos {
mat `var2' = nullmat(`var2') \ `var1'[`cn',1...]
}
mat colnames `wbeta' = `wendo'
mat `var_wbeta'=`var2' // and var_wbeta is the corresponding VCV
_estimates hold `waldmodel'
****************************** GET OPTION SPECS *************************************
* Now all varlists and counts are complete, so get option specs
get_option_specs, ///
`opts' ///
model(`model') /// pass info on model along with options
iid(`iid') /// pass iid boolean along with options
overid(`overid') /// pass overid boolean along with options
nendog(`nendog') /// pass #endog info along with options
nwendog(`nwendog') ///
nsendog(`nsendog') /// #strongly IDed endog; also boolean
ncsendog(`ncsendog') /// #not in subset test
npwendog(`npwendog') /// #weakd IDed for projection-based inference; also boolean
wbeta(`wbeta') /// row vector
var_wbeta(`var_wbeta'))
local level "`r(level)'"
local levellist "`r(levellist)'"
local ar_level "`r(ar_level)'"
local wald_level "`r(wald_level)'"
local clr_level "`r(clr_level)'"
local k_level "`r(k_level)'"
local j_level "`r(j_level)'"
local kj_level "`r(kj_level)'"
local kjk_level "`r(kjk_level)'"
local kjj_level "`r(kjj_level)'"
local clrsims "`r(clrsims)'"
local testlist "`r(testlist)'" // testlist for table
local citestlist "`r(citestlist)'" // testlist for CI table
local ptestlist "`r(ptestlist)'" // testlist for projection-based inference
local gridcols "`r(gridcols)'"
local gridmult "`r(gridmult)'"
local gridmin "`r(gridmin)'"
local gridmax "`r(gridmax)'"
local gridpoints "`r(gridpoints)'"
local null "`r(null)'"
local wnull "`r(wnull)'"
local kwt "`r(kwt)'"
local lm "`r(lm)'" // =0 or 1
local cuepoint "`r(cuepoint)'" // =0 or 1
local cuestrong "`r(cuestrong)'" // =0 or 1
local ci "`r(ci)'" // =0 or 1
local usegrid "`r(usegrid)'" // =0 or 1
local closedform "`r(closedform)'" // =0 or 1
local testid "`r(testid)'" // =0 or 1
local gridlist "`r(gridlist)'"
local gridpoints "`r(gridpoints)'"
local exportmats "`r(exportmats)'"
local eststorewald "`r(eststorewald)'"
local estadd =r(estadd) // =0 or 1
local estaddname "`r(estaddname)'"
local graph "`r(graph)'"
local graphxrange "`r(graphxrange)'"
local graphopt "`r(graphopt)'"
local contouropt "`r(contouropt)'"
local surfaceopt "`r(surfaceopt)'"
local contouronly "`r(contouronly)'"
local surfaceonly "`r(surfaceonly)'"
local displaywald "`r(displaywald)'"
local forcerobust "`r(forcerobust)'" // boolean
// end assembly of model and option specs
***************************** XT TRANSFORMS *************************************************
* Transform data (FE or FD) if required. Data already preserved.
if "`xtmodel'" == "fe" {
capture xtset
if _rc > 0 {
di as err "Fixed-effects estimation requires data to be -xtset-"
exit 459
}
local ivar "`r(panelvar)'"
local tvar "`r(timevar)'"
tempvar T_i
sort `ivar' `touse'
* Catch singletons. Must use unweighted data
qui by `ivar' `touse': gen long `T_i' = _N if _n==_N & `touse'
qui replace `touse'=0 if `T_i'==1
drop `T_i'
qui {
sort `ivar' `touse'
* Only iw and fw use weighted observation counts
if "`weight'" == "iweight" | "`weight'" == "fweight" {
by `ivar' `touse': gen long `T_i' = sum(`wvar') if `touse'
}
else {
by `ivar' `touse': gen long `T_i' = _N if `touse'
}
by `ivar' `touse': replace `T_i' = . if _n~=_N
* Demean. Create new vars as doubles (reusing existing non-double vars loses accuracy)
* Apply to TS temp vars if any so use _t list
local allvars_t "`depvar_t' `inexog_t' `endo_t' `exexog_t'"
foreach var of local allvars_t {
tempvar `var'_m
* To get weighted means
by `ivar' `touse' : gen double ``var'_m'=sum(`var'*`wvar')/sum(`wvar') if `touse'
by `ivar' `touse' : replace ``var'_m'=``var'_m'[_N] if `touse' & _n<_N
* This guarantees that the demeaned variables are doubles and have correct names
by `ivar' `touse' : replace ``var'_m'=`var'-``var'_m'[_N] if `touse'
drop `var'
rename ``var'_m' `var'
}
}
* Restore xtset-ing of data before leaving transform block
qui xtset `ivar' `tvar'
}
if "`xtmodel'" == "fd" {
* No transformation required!
* Handled by fvrevar above - varlists have D. operator in them.
* Code below is in case ivar or tvar are ever required.
capture xtset
if _rc > 0 {
di as err "Fixed-effects estimation requires data to be -xtset-"
exit 459
}
local ivar "`r(panelvar)'"
local tvar "`r(timevar)'"
}
******************************* PARTIAL-OUT *************************************************
* Always partial out exogenous regressors in linear models.
* `npartial' is #vars partialled out; includes constant in count
* `partialcons' =0 if no constant partialled out, =1 if it was
* `ninexog' = 0 if inexog all partialled out, otherwise unchanged
* `cons' = 0 if inexog all partialled out, otherwise unchanged
* `noconstant' = "noconstant" if exexog all partialled out, otherwise unchanged
* `inexog' = original varlist; unchanged
* `inexog_t' = empty if all inexog partialled out, otherwise unchanged
local npartial=0 // Default
local partialcons=0 // Default
if "`model'" == "linear" {
* Loop through variables to transform
* Use tempvars if any so that they are transformed
tempname partial_resid
foreach var of varlist `depvar_t' `endo_t' `exexog_t' {
qui regress `var' `inexog_t' if `touse' `wtexp', `noconstant' // `noconstant' here is from original model spec
qui predict double `partial_resid' if `touse', resid
qui replace `var' = `partial_resid'
drop `partial_resid'
}
* constant always partialled out
local noconstant "noconstant"
local partialcons =`cons'
local cons =0
* update count of included exogenous
local npartial =`ninexog'
local ninexog =0
* inexog_t list now empty but original inexog list unchanged
local inexog_t ""
} // end partial-out block
***************************** MISC PREPARATION ********************
* Small-sample adjustsment. Includes #partialled-out vars.
if "`small'"~="" & "`cluster'"=="" {
local ssa=(`N'-`dofminus')/(`N'-(`nexexog'+`ninexog'+`npartial')-`dofminus') // iid, robust, ac or HAC
}
else if "`small'"~="" & "`cluster'"~="" {
local ssa=(`N_clust')/(`N_clust'-1) * /// cluster
(`N'-1)/(`N'-(`nexexog'+`ninexog'+`npartial')) // no dofminus here
}
else {
local ssa=1 // no small
}
* Shared tempnames
tempname rk ar_p ar_chi2 ar_df k_p k_chi2 k_df j_p j_chi2 j_df kj_p ///
clr_p clr_stat clr_df ar_r k_r j_r kj_dnr kj_r clr_r ///
wald_p wald_chi2 wald_df wald_r
tempname nullvec del_z pi_z var_del var_pidel_z var_pi_z bhat del_v
tempname S S11 S12 S22 zz zzinv x2z xx zx xy x1x1 x2x2 x1x2 zx1 zx2 zy x1y x2y yy
tempname sbeta iv_sbeta0 pi2hat
tempname bhat pi_z1 pi_z2 var_pidel_z1 var_pidel_z2 var_pi_z11 var_pi_z12 var_pi_z22
tempname syy see sxx sxy sxe syv sve svv AA
tempname citable
tempname nullvector wnullvector
tempvar vhat vhat1 vhat2 uhat uhat1 uhat2 ehat
* Misc
local npd=0
******************************** CUE point estimator if requested **************************
if `cuepoint' {
tempname ecuebeta wcuebeta
if ~`overid' { // if exactly-ID, CUE=LIML=IV
mat `ecuebeta' = `ebeta'
mat `wcuebeta' = `wbeta'
}
else {
computecrossprods, touse(`touse') wtexp(`wtexp') exexog(`exexog_t') wendo(`endo_t') depvar(`depvar_t')
mat `zz' = r(zz)
mat `xx' = r(x1x1)
mat `zx' = r(zx1)
mat `xy' = r(x1y)
mat `zy' = r(zy)
mat `yy' = r(yy)
if `iid' {
mata: s_liml( ///
`N', ///
"`zz'", ///
"`xx'", ///
"`zx'", ///
"`xy'", ///
"`zy'", ///
"`yy'" ///
)
}
else {
tempvar cuewvar ehat
qui gen double `ehat' = . if `touse'
local avarcmd "avar (`ehat') (`exexog_t') if `touse' `wtexp', `vceopt' nocons"
qui gen double `cuewvar' = `wf'*`wvar' if `touse'
di as text "Obtaining CUE point estimates..."
mata: s_cue_beta( ///
`N', ///
"`depvar_t'", ///
"`endo_t'", ///
"`exexog_t'", ///
"`ehat'", ///
"`touse'", ///
"`ebeta'", /// use all Wald endog for starting values
"`cuewvar'", ///
"`zz'", ///
"`avarcmd'", ///
`lm', ///
"on" /// show CUE trace log
)
}
mat `ecuebeta'=r(beta) // CUE estimate for all endogenous
mat `ecuebeta'=`ecuebeta'' // row vector (Stata convention)
mat colnames `ecuebeta' = `endo'
foreach vn of local wendo { // CUE estimate for weakly-ID only
mat `wcuebeta' = nullmat(`wcuebeta') , `ecuebeta'[1,"`vn'"]
}
} // end CUE (overid) block
} // end cuepoint block
***************** PREPARE VECTOR OF NULLS INCLUDING PREP FOR WEAK/STRONG ************************
if "`wnull'"=="" {
mat `wnullvector' = J(1,`nwendog',0)
forvalues i=1/`nwendog' {
local wnull "`wnull' 0"
}
local wnull : list clean wnull
}
else {
foreach num in `wnull' {
mat `wnullvector' = nullmat(`wnullvector') , `num'
}
}
mat `nullvector' = `wnullvector' // null vector is weak nullvector
// followed by coeffs for strongly IDed (added later)
* Check that #nulls matches #weak endog
if `nwendog' ~= colsof(`wnullvector') {
di as err "error: number of weakly endogenous (`nwendog') doesn't match number of hypothesized values (" colsof(`wnullvector') ")"
exit 198
}
* weak/strong estimation - need efficient estimator of strongly-identified coeffs under null.
* Currently implemented for linear models only so inexog and constant already partialled out.
if `nsendog' > 0 {
* x1 is weakly identified (wendo), x2 is strongly identified (sendo)
computecrossprods, touse(`touse') wtexp(`wtexp') exexog(`exexog_t') wendo(`wendo_t') other(`sendo_t') depvar(`depvar_t')
mat `zz' = r(zz)
mat `x1x2' = r(x1x2)
mat `x1x1' = r(x1x1)
mat `x2x2' = r(x2x2)
mat `zx1' = r(zx1)
mat `zx2' = r(zx2)
mat `x1y' = r(x1y)
mat `x2y' = r(x2y)
mat `zy' = r(zy)
mat `yy' = r(yy)
if `iid' & `cuestrong' {
local s1method "liml"
}
else if `iid' {
local s1method "iv"
}
else if `cuestrong' {
local s1method "iv"
local s2method "cue"
}
else {
local s1method "iv"
local s2method "gmm2s"
}
* Always need to use an iid method (LIML or IV), either for final strong beta or for initial value for non-iid strong beta
get_strong_beta, ///
`s1method' /// either LIML or IV
nobs(`N') ///
b0(`wnullvector') ///
zz(`zz') ///
x1x1(`x1x1') ///
x1x2(`x1x2') ///
x2x2(`x2x2') ///
zx1(`zx1') ///
zx2(`zx2') ///
x1y(`x1y') ///
x2y(`x2y') ///
zy(`zy') ///
yy(`yy')
mat `sbeta' = r(sbeta)
* Next 2 will be missing if LIML
mat `iv_sbeta0' = r(iv_sbeta0) // IV beta for strongly-idenfified send at null=0 (used in grid search)
mat `pi2hat' = r(pi2hat) // RF coeffs matrix (used in grid search)
if ~`iid' { //
get_strong_beta, ///
`s2method' /// gmm2s or cue
lm(`lm') ///
nobs(`N') ///
b0(`wnullvector') ///
sbeta(`sbeta') /// 1st-step IV estimator for strong endog at specified null
zz(`zz') ///
x1x2(`x1x2') ///
zx1(`zx1') ///
zx2(`zx2') ///
zy(`zy') ///
depvar(`depvar_t') ///
wendo(`wendo_t') ///
sendo(`sendo_t') ///
exexog(`exexog_t') ///
touse(`touse') ///
vceopt(`vceopt') ///
wtexp(`wtexp') ///
wvar(`wvar') ///
wf(`wf')
mat `sbeta' = r(sbeta) // new sbeta based on efficient GMM
}
mat colnames `sbeta' = `sendo'
mat `nullvector' = `wnullvector' , `sbeta' // append coeffs for strong to nullvector
}
************************************ TEST SPECIFIED NULL *****************************
* Wald test. Works for all cases.
mata: s_wald( ///
"`var_wbeta'", ///
"`wbeta'", /// row vector
"`wnullvector'" /// rowvector
)
local wald_chi2 =r(wald_chi2)
local npd =max(`npd',r(npd)) // promote to 1 if npd matrix ever encountered
if `ncsendog' { /// subset AR test; requires IID and linearity
*** inexog and cons partialled out everywhere ***
* x1 is subset weakly identified (wendo), x2 complement (not in subset, not tested, coeffs obtained by LIML)
computecrossprods, touse(`touse') wtexp(`wtexp') exexog(`exexog_t') wendo(`wendo_t') other(`csendo_t') depvar(`depvar_t')
mat `zz' = r(zz)
mat `x1x2' = r(x1x2)
mat `x1x1' = r(x1x1)
mat `x2x2' = r(x2x2)
mat `zx1' = r(zx1)
mat `zx2' = r(zx2)
mat `x1y' = r(x1y)
mat `x2y' = r(x2y)
mat `zy' = r(zy)
mat `yy' = r(yy)
* subset AR is simply the overid stat (LM=Basmann, Wald/MD=Sargan) from LIML estimation
* with y0 = y - x1*nullvector, x1=weak endog, x2=other endog
mata: s_sliml( ///
`N', ///
"`zz'", ///
"`x1x1'", ///
"`x1x2'", ///
"`x2x2'", ///
"`zx1'", ///
"`zx2'", ///
"`x1y'", ///
"`x2y'", ///
"`zy'", ///
"`yy'", ///
"`wnullvector'", /// rowvector
`nexog', /// L = #exexog + #inexog; needed for AR but not LIML beta
0, /// flag=1 => calculate beta, =0 => calc only lambda and ar stat
`lm' /// flag LM=0, MD=1; needed for AR but not LIML beta
)
} // end subset-AR code
else if `iid' & ~`forcerobust' { // iid-based test formulae
computematrices_iid , ///
cons(`cons') /// cons=1 if constant in model (ivprobit/tobit)
touse(`touse') ///
wtexp(`wtexp') ///
model(`model') ///
depvar(`depvar_t') ///
endo(`wendo_t' `sendo_t') /// first weak, then strong
exexog(`exexog_t') ///
inexog(`inexog_t') /// empty/partialled-out unless ivprobit/tobit
nendog(`nendog') ///
nexexog(`nexexog') ///
ntinexog(`ntinexog') /// #exogenous regressors included in tests
npartial(`npartial') /// #inexog partialled-out
nobs(`N') ///
dofminus(`dofminus') ///
ssa(`ssa') /// small-sample adjustment
lm(`lm') ///
llopt(`llopt') /// ivtobit options
ulopt(`ulopt') ///
asis(`asis') // ivprobit options
mat `var_pi_z' =r(var_pi_z)
mat `pi_z' =r(pi_z)
mat `var_del' =r(var_del)
mat `del_z' =r(del_z)
mat `del_v' =r(del_v)
mat `syy' =r(syy)
mat `see' =r(see)
mat `sxy' =r(sxy)
mat `sve' =r(sve)
mat `sxx' =r(sxx)
mat `svv' =r(svv)
* Test specified null
mata: compute_tests( ///
1, /// 1=iid code, 0=robust code
`lm', ///
`nendog', ///
`nexexog', ///
"`del_z'", /// row vector
"`var_del'", ///
"`pi_z'", ///
"`var_pi_z'", ///
"`var_pidel_z'", ///
"`del_v'", /// row vector
"`wnullvector'", /// row vector
"`nullvector'", /// row vector; first weak, then strong (if any)
"`syy'", /// actually scalar
"`see'", /// actually scalar
"`sxy'", /// COLUMN VECTOR
"`sve'", /// COLUMN VECTOR
"`sxx'", /// KxK matrix
"`svv'" /// KxK matrix
)
} // end iid code
else { // robust/non-iid case
*** inexog and cons partialled out everywhere ***
computematrices_robust, ///
touse(`touse') ///
wtexp(`wtexp') ///
vceopt(`vceopt') ///
depvar(`depvar_t') ///
endo(`wendo_t' `sendo_t') /// first weak, then strong (if any)
exexog(`exexog_t') ///
nendog(`nendog') ///
nexexog(`nexexog') ///
npartial(`npartial') /// #inexog partialled-out
nobs(`N') ///
dofminus(`dofminus') ///
ssa(`ssa') /// small-sample adjustment
lm(`lm')
mat `del_z' = r(del_z)
mat `var_pi_z' = r(var_pi_z)
mat `pi_z' = r(pi_z)
mat `var_del' = r(var_del)
mat `var_pidel_z' = r(var_pidel_z)
* Test specified null
mata: compute_tests( ///
0, /// 1=iid code, 0=robust code
`lm', ///
`nendog', ///
`nexexog', ///
"`del_z'", /// row vector
"`var_del'", ///
"`pi_z'", ///
"`var_pi_z'", ///
"`var_pidel_z'", ///
"`del_v'", /// row vector
"`wnullvector'", /// row vector
"`nullvector'" /// row vector; first weak, then strong (if any)
)
} // end robust/non-iid code
************************** ASSEMBLE RESULTS FOR SPECIFIED NULL *****************************
* Collect test stats left behind by above
* If not overid or if subset, most will be missing
local ar_chi2 =r(ar_chi2)
local k_chi2 =r(k_chi2)
local j_chi2 =r(j_chi2)
local clr_stat =r(clr_stat)
local rk =r(rk)
local npd =max(`npd',r(npd)) // promote to 1 if npd matrix ever encountered
* calculate test statistics, p-values, and rejection indicators from above matrices
compute_pvals, ///
closedform(`closedform') ///
clrsims(`clrsims') ///
overid(`overid') ///
iid(`iid') ///
rk(`rk') ///
nexexog(`nexexog') ///
nendog(`nendog') ///
nsendog(`nsendog') ///
ncsendog(`ncsendog') ///
ar_level(`ar_level') ///
k_level(`k_level') ///
j_level(`j_level') ///
kj_level(`kj_level') ///
clr_level(`clr_level') ///
wald_level(`wald_level') ///
ar_chi2(`ar_chi2') ///
k_chi2( `k_chi2' ) ///
j_chi2(`j_chi2') ///
clr_stat(`clr_stat') ///
wald_chi2(`wald_chi2') ///
kwt(`kwt')
local wald_df =r(wald_df)
local wald_p =r(wald_p)
local wald_r =r(wald_r)
local ar_df =r(ar_df)
local ar_p =r(ar_p)
local ar_r =r(ar_r)
local clr_p =r(clr_p)
local clr_r =r(clr_r)
local kj_p =r(kj_p)
local kj_r =r(kj_r)
local k_df =r(k_df)
local k_p =r(k_p)
local k_r =r(k_r)
local j_df =r(j_df)
local j_p =r(j_p)
local j_r =r(j_r)
local rk_p =r(rk_p)
local rk_df =r(rk_df)
******************** Anderson/Kleibergen-Paap LM test of underid of Wald model *************
* Need ranktest installed
if `testid' {
capture ranktest, version
local vernum "`r(version)'"
if (_rc == 0) & ("`vernum'" >= "`ranktestversion'") {