/* 
HOUSEHOLD INCOMES IN TAX DATA: USING ADDRESSES TO MOVE FROM TAX UNIT TO HOUSEHOLD INCOME DISTRIBUTIONS
by Jeff Larrimore, Jacob Mortenson, and David Splinter

Data from IRS administrative data. 
*/ 


%macro grp995(dsn, var, wts, incgrpvar, p995_cutoff);
 * grp995 divides the sample into two, above and below the 99.5% level;
                    /* calculate the cutpoint */
                    proc univariate noprint data = &dsn;
                      var &var;
                      output out = income_grp pctlpts = 99.5 pctlpre = pct;
                      weight &wts; run;

                    /* write the quintiles to macro variables */
                    data _null_;  set income_grp;
                    call symput('cut',pct99_5); run;

                    /* create the new variable in the main dataset */
                    data &dsn;  set &dsn;
                      if &var =. then &incgrpvar = .;
                      if &var le &cut then &incgrpvar = .;
                      else &incgrpvar = &var;
                      &p995_cutoff = &cut; run;
%mend grp995;


%macro grp999(dsn, var, wts, incgrpvar, p999_cutoff);
 * grp999 divides the sample into two, above and below the 99.9% level;
                    /* calculate the cutpoint */
                    proc univariate noprint data = &dsn;
                      var &var;
                      output out = income_grp pctlpts = 99.9 pctlpre = pct;
                      weight &wts; run;

                    /* write the quintiles to macro variables */
                    data _null_;  set income_grp;
                    call symput('cut',pct99_9); run;

                    /* create the new variable in the main dataset */
                    data &dsn;  set &dsn;
                      if &var =. then &incgrpvar = .;
                      if &var le &cut then &incgrpvar = .;
                      else &incgrpvar = &var;
                      &p999_cutoff = &cut; run;
%mend grp999;


%macro centiles(dsn, var, wts, cent);
                    /* calculate the cutpoints for quintiles */
                    proc univariate noprint data = &dsn; var &var;
                      output out = centile pctlpts=1 to 99  pctlpre = pct; 
                      weight &wts; run;

                    /* write the decile cutoffs to macro variables */
                    data _null_; set centile;
*call symput('d0' , pct0) ;
call symput('d1' , pct1) ;
call symput('d2' , pct2) ;
call symput('d3' , pct3) ;
call symput('d4' , pct4) ;
call symput('d5' , pct5) ;
call symput('d6' , pct6) ;
call symput('d7' , pct7) ;
call symput('d8' , pct8) ;
call symput('d9' , pct9) ;
call symput('d10' , pct10);
call symput('d11' , pct11);
call symput('d12' , pct12);
call symput('d13' , pct13);
call symput('d14' , pct14);
call symput('d15' , pct15);
call symput('d16' , pct16);
call symput('d17' , pct17);
call symput('d18' , pct18);
call symput('d19' , pct19);
call symput('d20' , pct20);
call symput('d21' , pct21);
call symput('d22' , pct22);
call symput('d23' , pct23);
call symput('d24' , pct24);
call symput('d25' , pct25);
call symput('d26' , pct26);
call symput('d27' , pct27);
call symput('d28' , pct28);
call symput('d29' , pct29);
call symput('d30' , pct30);
call symput('d31' , pct31);
call symput('d32' , pct32);
call symput('d33' , pct33);
call symput('d34' , pct34);
call symput('d35' , pct35);
call symput('d36' , pct36);
call symput('d37' , pct37);
call symput('d38' , pct38);
call symput('d39' , pct39);
call symput('d40' , pct40);
call symput('d41' , pct41);
call symput('d42' , pct42);
call symput('d43' , pct43);
call symput('d44' , pct44);
call symput('d45' , pct45);
call symput('d46' , pct46);
call symput('d47' , pct47);
call symput('d48' , pct48);
call symput('d49' , pct49);
call symput('d50' , pct50);
call symput('d51' , pct51);
call symput('d52' , pct52);
call symput('d53' , pct53);
call symput('d54' , pct54);
call symput('d55' , pct55);
call symput('d56' , pct56);
call symput('d57' , pct57);
call symput('d58' , pct58);
call symput('d59' , pct59);
call symput('d60' , pct60);
call symput('d61' , pct61);
call symput('d62' , pct62);
call symput('d63' , pct63);
call symput('d64' , pct64);
call symput('d65' , pct65);
call symput('d66' , pct66);
call symput('d67' , pct67);
call symput('d68' , pct68);
call symput('d69' , pct69);
call symput('d70' , pct70);
call symput('d71' , pct71);
call symput('d72' , pct72);
call symput('d73' , pct73);
call symput('d74' , pct74);
call symput('d75' , pct75);
call symput('d76' , pct76);
call symput('d77' , pct77);
call symput('d78' , pct78);
call symput('d79' , pct79);
call symput('d80' , pct80);
call symput('d81' , pct81);
call symput('d82' , pct82);
call symput('d83' , pct83);
call symput('d84' , pct84);
call symput('d85' , pct85);
call symput('d86' , pct86);
call symput('d87' , pct87);
call symput('d88' , pct88);
call symput('d89' , pct89);
call symput('d90' , pct90);
call symput('d91' , pct91);
call symput('d92' , pct92);
call symput('d93' , pct93);
call symput('d94' , pct94);
call symput('d95' , pct95);
call symput('d96' , pct96);
call symput('d97' , pct97);
call symput('d98' , pct98);
call symput('d99' , pct99);
 run; 
					* use centile cutoffs to set centile variable cent;
                    data &dsn; set &dsn;
                      if &var =. then &cent =.;
					  %do centile = 1 %to 99; 
                            if (&cent=.) and (&var < &&d&centile) then &cent=&centile; 
					  %end;
					   if (&var > &d99) then &cent=100;
                      run;
%mend centiles;


%macro gini3(capinc, capwts, dsin, yr, var); 
/*  This SAS code was originally written by Philip N. Cohen.  */

title 'Income distribution';
proc freq data = &dsin;
tables &capinc / noprint out = table;
format &capinc 7.0;
weight &capwts;
run;

data table; set table;
retain suminc perpop;

suminc + (&capinc * count);
perpop + percent;

run;


proc sort data=table;
by descending suminc; run;

data table; set table;
by descending suminc;

if _n_=1 then do;
totalinc=suminc; end;

retain totalinc;
perinc = (suminc/totalinc) * 100;
run;

proc sort data=table;
by perpop; run;

/* To calculate Gini:
   sum[Xsub(i) * Ysub(i+1)] - sum[Xsub(i+1) * Ysub(i)]
   where X is the proportion of population column and
   Y is the proportion of income column.  */
data ginidat; set table;

xlag = lag(perpop);
xlag = xlag / 100;

ylag = lag(perinc);
ylag = ylag / 100;

columna = (perinc/100) * xlag;
columnb = (perpop/100) * ylag;

retain suma sumb;

suma + columna;
sumb + columnb;
gini&yr = suma - sumb;

year = &yr;
gini&var = gini&yr;
if (perinc ne 100) then delete;
keep year gini&var;
run;

data gini&var; set ginidat; keep year gini&var; run;

%mend gini3;


%macro TaxUnitDist2010;

%let NF2010 = 9663;

data nf_txunits; set INDIV.nonfiler10; * JCT ITM extract of nonfilers;
  gross_ret = ptpen + tirad; * ITM splits 1099-R income into pension and IRA amounts, tirad is gross IRA distributions minus rollovers;
  roll = 0;
  wages = was;
  div = dbe;
  intr = intst; * missing nontaxable interest;
  misc = tot_misc_inc;
  ssa = ssinc;
  ui = ucagi; * this is all UI on 1099-G;
  xwgt=wt;

  HHsize=1+ depne; if (mars=2) then HHsize=HHsize+1; * if married then add one more person to tax unit (about 40% nonfiling tax units are married);
  if (HHsize=.)or(HHsize<1) then HHsize=1; 
  HHwt = xwgt*HHsize;

  incmkt  = wages + div + intr + 0.3*misc; * Market income (PS comparison);
  inc5adj = wages + div + intr + 0.3*misc + ssa + ui; * Basic income (Census comparison);
  incnew =  wages + div + intr + 0.3*misc + ssa + ui + gross_ret - roll; * Pre-tax income: non-filer tax units;
  inctx  =  wages + div + intr + 0.3*misc + ssa + ui + gross_ret - roll; * Private income plus SS/DI/UI/TAX cash transfers = Private income plus SS/DI/UI cash transfers + Refundable EIC/CTC;
  incat  = inctx ; * assume non-filers pay no federal tax and get no refundable credits;

  szadj_new   = incnew;  if (HHsize>0) then szadj_new   = incnew/sqrt(HHsize);  if (szadj_new<1)  and(szadj_new ne .)   then szadj_new  =1+ranuni(364); 
  szadj_i5adj = inc5adj; if (HHsize>0) then szadj_i5adj = inc5adj/sqrt(HHsize); if (szadj_i5adj<1)and(szadj_i5adj ne .) then szadj_i5adj=1+ranuni(364); 
  szadj_tx    = inctx; if (HHsize>0) then szadj_tx = inctx/sqrt(HHsize);        if (szadj_tx<1)   and(szadj_tx ne .)    then szadj_tx   =1+ranuni(364);
  szadj_at    = incat; if (HHsize>0) then szadj_at = incat/sqrt(HHsize);        if (szadj_at<1)   and(szadj_at ne .)    then szadj_at   =1+ranuni(364);

  if (incmkt <1)and(incmkt  ne .) then incmkt =1+ranuni(364);
  if (inc5adj<1)and(inc5adj ne .) then inc5adj=1+ranuni(364);
  if (incnew <1)and(incnew  ne .) then incnew =1+ranuni(364); 
  if (inctx  <1)and(inctx   ne .) then inctx  =1+ranuni(364); 
  if (incat  <1)and(incat   ne .) then incat  =1+ranuni(364); 

  ct = 1;
  if (misc<0) then delete; * remove nonfiling tax units with foreign income exclusion, note that these were based on 2007 filers, but with 2010 data and this was the reconciliation method;

  keep incnew incmk inc5adj inctx incat szadj_i5adj szadj_tx szadj_at xwgt ct mars HHsize HHwt szadj_new;
run;

* Filing Tax Units for 2010 (actually these are tax returns filed in calendar year 2011, and so includes a few percent late filers);
data taxunits; set indiv.soi2010;
   id         = S002; * primary filer TIN;
   ids        = S003; * spouses TIN;   
   flpdyr = flpdyr;  
   EXEM = N3+N6;           IF EXEM = . THEN EXEM = 1;       * Taxpayer+Dependents Exemptions;
   NotLive = N8;           IF NotLive = . THEN NotLive = 0; * Taxpayer+Dependents Exemptions;
   INCM = E00050;          IF INCME = . THEN INCME = 0;
   AGIX = E00100;          IF AGIX = . THEN AGIX = 0.0;
   wage_amt = E00200;      IF WAGE = . THEN WAGE = 0;
   INTEREST = E00300;      IF INTEREST = . THEN INTEREST = 0; * INTEREST RECEIVED, E253809 other interest (sch E), E00400 Tax-exempt interest, E21090 interest sch B;
   taxexint = E00400;      IF taxexint = . THEN taxexint = 0; * tax-exempt interest;
   DIVID = E00600;         IF DIVID = . THEN DIVID = 0;       * DIVIDEND INCOME (sch B);    
   cpftloss = E00900;      IF cpftloss = . THEN cpftloss = 0.0;  * BUSINESS/PROF NET PROFIT/LOSS + Combined Partnership and S-Corp net income or loss + SCH E NET INCOME OR LOSS ?????????????;
   PARTSCP = E26270;       IF PARTSCP = . THEN PARTSCP = 0;
   dpftloss = E01000;      IF dpftloss = . THEN dpftloss = 0.0;  * NET CAP GN/LOSS SCH D; 
   CGDS = E01100;          IF CGDS = . THEN CGDS = 0.0;          * E01100 CAP GN DISTRIB (1040) + E01200 supplemental sch net gain/loss;
   PENSION = E01300+E01500;           IF PENSION = . THEN PENSION = 0;  * E01500 Total pensions and annuities received(Form) (E01700 Pensions and annuities in AGI);
   PENSNT =  (E01300-E01400)+E01800;  IF PENSNT = . THEN PENSNT = 0;    * Nontaxable pensions and annuities (not in AGI);
   epftloss = E02000;      IF epftloss=. THEN epftloss=0; * Sch E;
   fpftloss = E02100;      IF fpftloss = . THEN fpftloss  = 0.0;   * FARM NET PROFIT OR LOSS + E58500 Net farm rental income (loss) from 4835;
   UNEM = E02300;          IF UNEM = . THEN UNEM = 0.0;            * Unemp compensation in AGI;
   txsocsec = E02500;      IF txsocsec = . THEN txsocsec = 0.0; 
   ssa = E02400;           IF ssa = . THEN ssa = 0.0; * taxable and non-taxable SS income;
   tot_inc = E02650;       IF tot_inc = . THEN tot_inc = 0.0; 
   STATADJ = E02900;       IF STATADJ = . THEN STATADJ = 0.0; * Total statutory adjustments (moving, payments to IRA, Keogh, alimony, etc.);
   ACTC   = e11070;        IF (ACTC = .)or(ACTC<0) THEN ACTC = 0.0;* e11070 | Additional Child Tax Credit; 
   TAXES  = e10300;   * line 60, Total tax liability ;
   EIC    = e59660;        IF (EIC = .)or(EIC<0) THEN EIC = 0.0; * e59660 (form value is e11000) | Total earned income credit; 
   taxIRA = E01400;        IF taxIRA = . THEN taxIRA = 0.0;      * taxable IRA distributions;
   taxPEN = E01700;        IF taxPEN = . THEN taxPEN = 0;        * taxable pensions and annuities;
   txret = taxIRA+taxPEN;

  if (txsocsec > ssa) then ssa = txsocsec; * often leave ssa blank;
  HHsize = EXEM; * Tax Unit size is total nubmer of exemptions (includes dep children not living with filer as in other data);
  if (HHsize=.)or(HHsize<1) then HHsize=1; 
  HHwt = xwgt*HHsize;

  incmk   = tot_inc - dpftloss - txsocsec - UNEM; * PS market income;
  inc5adj = tot_inc + taxexint - dpftloss - txret - txsocsec + ssa; * Basic Income = Private income less retirement (keeping SS) = 1040 total income + tax exempt interest – Sch D cap gains – taxable retirement (IRA+pensions);  
  incnew  = tot_inc + taxexint - dpftloss         - txsocsec + ssa; * Pre-tax income: filers = 1040 total income + tax-exempt interest + non-taxable retirement (IRA distributions, pensions and annuities) less realized capital gains; 
  inctx   = tot_inc + taxexint - txsocsec + ssa ; * Private income plus SS/DI/UI/TAX cash transfers = Private income plus SS/DI/UI cash transfers + Refundable EIC/CTC;
  incat = inctx + EIC + ACTC - TAXES; * Federal after-tax Private income plus SS/DI/UI/TAX cash transfers = Private income plus SS/DI/UI cash transfers + Refundable EIC/CTC – Federal Income Tax;

  if (incmkt <1)and(incmkt  ne .) then incmkt =1+ranuni(364);
  if (inc5adj<1)and(inc5adj ne .) then inc5adj=1+ranuni(364);
  if (incnew <1)and(incnew  ne .) then incnew =1+ranuni(364); 
  if (inctx  <1)and(inctx   ne .) then inctx  =1+ranuni(364); 
  if (incat  <1)and(incat   ne .) then incat  =1+ranuni(364);  

  szadj_new   = incnew;  if (HHsize>0) then szadj_new   = incnew /sqrt(HHsize); 
  szadj_i5adj = inc5adj; if (HHsize>0) then szadj_i5adj = inc5adj/sqrt(HHsize); 
  szadj_tx    = inctx;   if (HHsize>0) then szadj_tx    = inctx  /sqrt(HHsize);                  
  szadj_at    = incat;   if (HHsize>0) then szadj_at    = incat  /sqrt(HHsize); 

  szadj_new   = incnew;  if (HHsize>0) then szadj_new   = incnew/sqrt(HHsize);  if (szadj_new<1)  and(szadj_new ne .)   then szadj_new  =1+ranuni(364); 
  szadj_i5adj = inc5adj; if (HHsize>0) then szadj_i5adj = inc5adj/sqrt(HHsize); if (szadj_i5adj<1)and(szadj_i5adj ne .) then szadj_i5adj=1+ranuni(364); 
  szadj_tx    = inctx;   if (HHsize>0) then szadj_tx = inctx/sqrt(HHsize);      if (szadj_tx<1)   and(szadj_tx ne .)    then szadj_tx   =1+ranuni(364);
  szadj_at    = incat;   if (HHsize>0) then szadj_at = incat/sqrt(HHsize);      if (szadj_at<1)   and(szadj_at ne .)    then szadj_at   =1+ranuni(364);

  incat_norc = inctx - TAXES; * After-tax but no ref. credits;
  if (incat_norc<1)and(incat_norc ne .) then incat_norc =1+ranuni(364);

  rc=.; if (EIC + ACTC)>0 then rc=EIC + ACTC;
  rc2=0; if (EIC + ACTC)>0 then rc2=EIC + ACTC;
  refcr = 0; if (rc ne .) then refcr=1;
  incat_rc =.; if (refcr=1) then incat_rc =incat;
  incat_rcz=.; if (refcr=1) then incat_rcz=incat_norc;

  eic_cr= 0; if (eic ne .)and(eic>0) then eic_cr=1;
  inctaxamt = inctx - incat; * inctaxamtadj = szadj_tx - szadj_at; 

  if (stim2007=1) then delete; * No dependents or 2007 stimulus filers for ITM blowups;
  if (dsi=1) then delete; * remove dependent filers as not really separate tax units, causes an underestimate of nonfilers;
  ct = 1;
  gp = 1; if (inc5adj>100000) and (inc5adj<350000) then gp=2;  if (inc5adj>350000) and (inc5adj<1000000) then gp=3;  if (inc5adj>1000000) then gp=4;
  keep id ids incnew incmk inc5adj inctx incat szadj_i5adj szadj_tx szadj_at xwgt gp ct mars HHsize HHwt rc rc2 refcr inctaxamt incat_rc incat_rcz incat_norc eic actc eic_cr szadj_new;
run;

 
* Figure 1: merge information return income to filers (no non-filers) and no size-adjustment here, see CWHS-HHshocks-OLD program for more general info return code;
libname IRMF2010 '...';

data w2; set IRMF2010.w2; 
  wagetip  = EW006; id = wc033; ids=id; if id>0;
  keep id ids wagetip; run;
proc sort; by id; run;
 * Aggregate by SSN;
data w2 (keep=id was ids); retain was; set w2;
	by id;
	if first.id then was=wagetip; * initialize wage;
	    else was=was+wagetip; * add wage from other employers;
	if last.id;
run;
* merge primary filers;
data w21; set w2; drop ids; run;
proc sort data=taxunits; by id; run;
proc sort data=w21;      by id; run;
data taxunits; merge taxunits (in=a) w21; by id; if a; run;
* merge secondary filers;
data w2; set w2; wass=was; keep ids wass; run;
proc sort data=taxunits; by ids; run;
proc sort data=w2;       by ids; run;
data taxunits; merge taxunits (in=a) w2; by ids; if a; run;

*1099-DIV (dividend income); 
 data div; set IRMF2010.f1099div; 
     id=wc033;  div=ew065; if div=. then div=0; ids=id;
     keep id ids div; run;
proc sort; by id; run;
 * Aggregate by SSN;
data div (keep=id ids divs); retain divs; set div;
	by id;
	if first.id then divs=div; * initialize dividend;
	    else divs=divs+div; * add dividends from others;
	if last.id;
run;
* merge primary filers;
data div1; set div; drop ids; run;
proc sort data=taxunits; by id; run;
proc sort data=div1;     by id; run;
data taxunits; merge taxunits (in=a) div1; by id; if a; run;
* merge secondary filers;
data div; set div; divss=divs; keep ids divss; run;
proc sort data=taxunits; by ids; run;
proc sort data=div;      by ids; run;
data taxunits; merge taxunits (in=a) div; by ids; if a; run;

*1099-G (certain government payments including unemployment compensation);
 data f1099g; set IRMF2010.f1099g; 
     id=wc033; if ew020=. then ew020=0; if ew001=. then ew001=0; gov=ew020-ew001; if gov=. then gov=0; ids=id;
     keep id ids gov; run;
proc sort; by id gov; run;
 * Aggregate by SSN;
data f1099g (keep=id ids govs); retain govs; set f1099g;
	by id;
	if first.id then govs=gov; * initialize ;
	    else govs=govs+gov;  * add from others;
	if last.id;
run;
data f1099g1; set f1099g; drop ids; run;
proc sort data=taxunits; by id; run;
proc sort data=f1099g1;  by id; run;
data taxunits; merge taxunits (in=a) f1099g1; by id; if a; run;
data f1099g; set f1099g; govss=govs; keep ids govss; run;
proc sort data=taxunits; by ids; run;
proc sort data=f1099g;   by ids; run;
data taxunits; merge taxunits (in=a) f1099g; by ids; if a; run;

*1099-INT (interest income);
 data int; set IRMF2010.f1099int;
     id=wc033; a1=ew001; a2=ew002; a3=ew034; if (a1=.) then a1=0; if (a2=.) then a2=0; if (a3=.) then a3=0; int=a1+a2+a3; ids=id;
     keep id ids int; run;
proc sort; by id int; run;
 * Aggregate by SSN and year;
data int (keep=id ids ints); retain ints; set int;
	by id;
	if first.id then ints=int; * initialize;
	    else ints=ints+int;  * add from others;
	if last.id;
run;
data int1; set int; drop ids; run;
proc sort data=taxunits; by id; run;
proc sort data=int1;     by id; run;
data taxunits; merge taxunits (in=a) int1; by id; if a; run;
data int; set int; intss=ints; keep ids intss; run;
proc sort data=taxunits; by ids; run;
proc sort data=int;      by ids; run;
data taxunits; merge taxunits (in=a) int; by ids; if a; run;

*1099-MISC (miscellaneous income);
 data misc; set IRMF2010.f1099misc;
     id=wc033; 
                               a21=ew021; a23=ew023; a25=ew025; a32=ew032; a48=ew048; a49=ew049; a61=ew061; a153=ew153; a172=ew172; a173=ew173; 
                               if (a21=.) then a21=0; if (a23=.) then a23=0; if (a25=.) then a25=0; if (a32=.) then a32=0; if (a48=.) then a48=0; if (a49=.) then a49=0; 
							   if (a61=.) then a61=0; if (a153=.) then a153=0; if (a172=.) then a172=0; if (a173=.) then a173=0; 
                               misc=a21+a23+a25+a32+a48+a49+a61+a153+a172+a173;  if misc=. then misc=0; ids=id;
     keep id ids misc; run;
proc sort; by id misc; run;
 * Aggregate by SSN and year;
data misc (keep=id ids miscs); retain miscs; set misc;
	by id;
	if first.id then miscs=misc; * initialize ;
	    else miscs=miscs+misc;  * add from others;
	if last.id;
run;
data misc1; set misc; drop ids; run;
proc sort data=taxunits; by id; run;
proc sort data=misc1;  by id; run;
data taxunits; merge taxunits (in=a) misc1; by id; if a; run;
data misc; set misc; miscss=miscs; keep ids miscss; run;
proc sort data=taxunits; by ids; run;
proc sort data=misc;     by ids; run;
data taxunits; merge taxunits (in=a) misc; by ids; if a; run;

*1099-R (retirement distributions);
 data ret; set IRMF2010.f1099r;
     id=wc033; ret=ew128; if ret=. then ret=0; ids=id; 
     keep id ids ret; run;
proc sort; by id ret; run;
 * Aggregate by SSN and year;
data ret (keep=id rets); retain rets; set ret;
	by id;
	if first.id then rets=ret; * initialize ;
	    else rets=rets+ret;  * add from others;
	if last.id;
run;
data ret1; set ret; drop ids; run;
proc sort data=taxunits; by id; run;
proc sort data=ret1;  by id; run;
data taxunits; merge taxunits (in=a) ret1; by id; if a; run;
data ret; set ret; retss=rets; keep ids retss; run;
proc sort data=taxunits; by ids; run;
proc sort data=ret;      by ids; run;
data taxunits; merge taxunits (in=a) ret; by ids; if a; run;

*1099-SSA (Social Security income);
 data f1099ssa; set IRMF2010.f1099ssa;
     id=wc033; a1=ew004; a2=ew120; if (a1=.) then a1=0; if(a2=.) then a2=0; ssai=a1+a2; if ssai=. then ssai=0; ids=id;
     keep id ids ssai; run;
proc sort; by id ssai; run;
 * Aggregate by SSN and year;
data f1099ssa (keep=id ids ssas); retain ssas; set f1099ssa;
	by id;
	if first.id then ssas=ssai; * initialize ;
	    else ssas=ssas+ssai;  * add from others;
	if last.id;
run;
data f1099ssa1; set f1099ssa; drop ids; run;
proc sort data=taxunits;  by id; run;
proc sort data=f1099ssa1; by id; run;
data taxunits; merge taxunits (in=a) f1099ssa1; by id; if a; run;
data f1099ssa; set f1099ssa; ssass=ssas; keep ids ssass; run;
proc sort data=taxunits; by ids; run;
proc sort data=f1099ssa; by ids; run;
data taxunits; merge taxunits (in=a) f1099ssa; by ids; if a; run;
 
* sum information return income and set misc to 30% to account for deductions;
data taxunits; set taxunits;
   if (was=.)   then was=0;   if (wass=.)   then wass=0;
   if (divs=.)  then divs=0;  if (divss=.)  then divss=0;
   if (govs=.)  then govs=0;  if (govss=.)  then govss=0;
   if (ints=.)  then ints=0;  if (intss=.)  then intss=0;
   if (miscs=.) then miscs=0; if (miscss=.) then miscss=0;
   if (rets=.)  then rets=0;  if (retss=.)  then retss=0;
   if (ssas=.)  then ssas=0;  if (ssass=.)  then ssass=0;
   inc_inforet = was + wass + divs + divss + govs + govss + ints + intss + 0.3*(miscs + miscss) + ssas + ssass;
   if (inc5adj<1) then inc5adj=1; if (inc_inforet<1) then inc_inforet=1; * bottom-code incomes at $1; 
   *if (inc_inforet>100000)and(incnew<25000) then inc_inforet=100000; * income limitation: if tax return income is under $25,000 then information return income is a maximum of $100,000;
run;

*Figure 1;
%centiles(taxunits, inc5adj, xwgt, cent80);
proc sort data=taxunits; by cent80; run; 
proc means noprint data=taxunits; by cent80; var inc5adj; weight xwgt;
         output out=TaxUnits_Fig1 sum(ct HHsize)=taxunits Indivs  mean(inc5adj inc_inforet)=meanninc_noret mean_inc_inforet; run;

proc means noprint data=taxunits; var xwgt; 
         output out=Filers_Tot  sum(xwgt)= weight; run;

proc sort data=taxunits; by mars gp; run; 
proc means noprint data=taxunits; by mars gp; var inc5adj; weight xwgt;
         output out=TaxUnits_byMars  max(inc5adj)= mxinc5 min(inc5adj)= mininc5 mean(inc5adj)= mninc5 sum(inc5adj)=sum; run;

proc means noprint data=nf_txunits; var xwgt; 
         output out=Nonfilers_Tot  sum(xwgt)= weight; run;

* Append nonfilers and filers;
proc append base=taxunits data=nf_txunits force; run; 

data taxunits; set taxunits; 
  if (rc2=.) then rc2=0; 
  if (refcr      = .) then refcr    =0; * nonfilers have no ref. credit;
  if (incat_rc   = .) then incat_rc =0;
  if (incat_norc = .) then incat_rcz=0; 
  if (incmk  <1)or(incmk   = .) then incmk  =1+ranuni(364); 
  if (inc5adj<1)or(inc5adj = .) then inc5adj=1+ranuni(364); 
  if (incnew <1)or(incnew  = .) then incnew =1+ranuni(364);
  if (inctx  <1)or(inctx   = .) then inctx  =1+ranuni(364); 
  if (incat  <1)or(incat   = .) then incat  =1+ranuni(364); 
  if (incmk   = &NF2010)        then incmk  =&NF2010+ranuni(364);
  if (inc5adj = &NF2010)        then inc5adj=&NF2010+ranuni(364);
  if (szadj_i5adj<1)and(szadj_i5adj ne .) then szadj_i5adj=1+ranuni(364); 
  if (szadj_tx<1)and(szadj_tx ne .) then szadj_tx=1+ranuni(364);   
  if (szadj_at<1)and(szadj_at ne .) then szadj_at=1+ranuni(364);
  if (szadj_new<1)and(szadj_new ne .) then szadj_new=1+ranuni(364);
run;

proc sort data=taxunits; by HHsize; run; 
proc means noprint data=taxunits; by HHsize; var ct; weight xwgt;
         output out=TaxUnits_SizeDist sum(ct)=ct; run;

* Size-adjusted Census-type income grouped by #individuals;
%centiles(taxunits, szadj_i5adj, HHwt, cent);
proc sort data=taxunits; by cent; run; 
proc means noprint data=taxunits; by cent; var szadj_i5adj; weight HHwt;
         output out=TaxUnits_bySzAdj5adj  sum(ct)=Indivs  min(szadj_i5adj)=min_inc  mean(szadj_i5adj)=mninc5  sum(szadj_i5adj)=sum_inc; run;


* Figures 5,6,7,8: SIZE-ADJUSTED INCOME;
%centiles(taxunits, szadj_new, HHwt, cent2);
proc sort data=taxunits; by cent2; run; 
proc means noprint data=taxunits; by cent2; var szadj_new; weight HHwt;
         output out=TaxUnits_SzAdj sum(ct HHsize)=taxunits Indivs  min(szadj_new)=min_inc mean(szadj_new)=mnszadj sum(szadj_new)=sum_inc  
                                   mean(eic_cr)=eic_cr; run;
proc means noprint data=taxunits; by cent2; var szadj_new; weight xwgt;
         output out=TaxUnits_EIC sum(ct HHsize)=taxunits Indivs  mean(eic)=eic; run;
 
* Figure 9: Average tax rates by pre-tax income (includes cap gains);
%centiles(taxunits, szadj_tx, HHwt, cent45);
proc sort data=taxunits; by cent45; run; 
proc means noprint data=taxunits; by cent45; var szadj_tx; weight xwgt; * tax unit weight so that only count credits and taxes once, CBO approach;
         output out=TaxRates_txunits sum(ct)=HHs mean(inctx incat szadj_tx szadj_at inctaxamt)=inctx incat szadj_tx szadj_at inctaxamt; run; 
	
* Table 3: Size-Adjusted Income;
%centiles(taxunits, szadj_new, HHwt, cent222);
%grp995(taxunits,   szadj_new, HHwt, p995, p995_cutoff);
%grp999(taxunits,   szadj_new, HHwt, p999, p999_cutoff);

data taxunits; set taxunits;
q1=.; q2=.; q3=.; q4=.; top8095=.; top5=.;
       if (cent222<=20) then q1=szadj_new; *incnew; 
  else if (cent222<=40) then q2=szadj_new; *incnew; 
  else if (cent222<=60) then q3=szadj_new; *incnew; 
  else if (cent222<=80) then q4=szadj_new; *incnew; 
  else if (cent222<=95) then top8095=szadj_new; *incnew; 
  else if (cent222 >95) then top5=szadj_new; *incnew;
  pos5=szadj_new; if (szadj_new<0) then pos5=0;
  q1pos=q1; if (incnew<0) then q1pos=0;
  top1=.; if (cent222=100) then top1=szadj_new; *incnew;
  year=2010;
run;

proc means noprint data=taxunits; var incnew; weight HHwt; * individual weights;
         output out=Ineq_TU_szadj_Shares mean(year)=year sum(ct szadj_new pos5 q1 q1pos q2 q3 q4 top8095 top5)=ct tot_incnew tot_pos q1 q1pos q2 q3 q4 top8095 top5
                             P10(szadj_new)=p10 p20(szadj_new)=p20 p40(szadj_new)=p40 p50(szadj_new)=p50 p60(szadj_new)=p60 p80(szadj_new)=p80 p90(szadj_new)=p90 p95(szadj_new)=p95 p99(szadj_new)=p99
                             mean(q1 q1pos q2 q3 q4  top8095 top5)=mn_q1 mn_q1pos mn_q2 mn_q3 mn_q4 mn_top8095 mn_top5 
                             sum(top1 p995 p999)=top1 top05 top01 min(p995 p999)=p995 p999 mean(top1 p995 p999)=mn_top1 mn_top05 mn_top01; run;

%gini3(szadj_new, HHwt, taxunits, 2010, _szadj); 

*Not size-adjusted;
*Basic income (census comparison) grouped by #tax units;
%centiles(taxunits, incnew, xwgt, cent2);
proc sort data=taxunits; by cent2; run; 
proc means noprint data=taxunits; by cent2; var incnew; weight xwgt;
         output out=TaxUnits_IncNew sum(ct HHsize)=taxunits Indivs  min(incnew)=min_inc  mean(incnew)=mninc5  sum(incnew)=sum_inc
            mean(eic_cr eic)=eic_cr eic; run;

* Average tax rates by pre-tax income (includes cap gains);
%centiles(taxunits, inctx, xwgt, cent30);
proc sort data=taxunits; by cent30; run; 
proc means noprint data=taxunits; by cent30; var inctx; weight xwgt;
         output out=TaxRates_txunits sum(ct)=HHs mean(inctx incat)=inctx incat; run; 

%centiles(taxunits, incnew, xwgt, cent222);
%grp995(taxunits,   incnew, xwgt, p995, p995_cutoff);
%grp999(taxunits,   incnew, xwgt, p999, p999_cutoff);

data taxunits; set taxunits;
q1=.; q2=.; q3=.; q4=.; top8095=.; top5=.;
       if (cent222<=20) then q1=incnew; 
  else if (cent222<=40) then q2=incnew; 
  else if (cent222<=60) then q3=incnew; 
  else if (cent222<=80) then q4=incnew; 
  else if (cent222<=95) then top8095=incnew; 
  else if (cent222 >95) then top5=incnew;
  pos5=incnew; if (incnew<0) then pos5=0;
  q1pos=q1; if (incnew<0) then q1pos=0;
  top1=.; if (cent222=100) then top1=incnew;
  year=2010;
run;

proc means noprint data=taxunits; var incnew; weight xwgt;
         output out=Ineq_TU_new mean(year)=year sum(ct incnew pos5 q1 q1pos q2 q3 q4 top8095 top5)=ct tot_incnew tot_pos q1 q1pos q2 q3 q4 top8095 top5
                             P10(incnew)=p10 p20(incnew)=p20 p40(incnew)=p40 p50(incnew)=p50 p60(incnew)=p60 p80(incnew)=p80 p90(incnew)=p90 p95(incnew)=p95 p99(incnew)=p99
                             mean(q1 q1pos q2 q3 q4  top8095 top5)=mn_q1 mn_q1pos mn_q2 mn_q3 mn_q4 mn_top8095 mn_top5 
                             sum(top1 p995 p999)=top1 top05 top01 min(p995 p999)=p995 p999 mean(top1 p995 p999)=mn_top1 mn_top05 mn_top01; run;

%gini3(incnew, xwgt, taxunits, 2010, _new);

%mend TaxUnitDist2010;


%macro CDWdata(yr);

 * 2010 file was called FuzzyMerge_2018_b or 180510;
proc import datafile="....extract_5pct.csv"
  out=hhs_in dbms=csv replace; guessingrows=50000; getnames=y; delimiter=','; * 6257177 observations;
run;
  
data hhs; set hhs_in; 
  if tot_inc=. then tot_inc=0; if dpftloss=. then dpftloss=0; if taxexint=. then taxexint=0; if txsocsec=. then txsocsec=0; if txret=. then txret=0; if wage_amt=. then wage_amt=0; if fpftloss=. then fpftloss=0;
  if cpftloss=. then cpftloss=0; if epftloss=. then epftloss=0; if ui=. then ui=0; if totximf=. then totximf=0; if eic=. then eic=0; if addlchc=. then addlchc=0;
  if ssa=. then ssa=0; if rollovers=. then rollovers=0; if gross_ret=. then gross_ret=0; if wages=. then wages=0; if div=. then div=0; if intr=. then intr=0; if misc=. then misc=0;
  if finf_wages=. then finf_wages=0; if finf_ssa=. then finf_ssa=0; if finf_un=. then finf_un=0; if finf_ui=. then finf_ui=0; if finf_misc=. then finf_misc=0; if pobox=. then pobox=0; 
  if multiyrtrunc=. then multiyrtrunc=0;  

  wt = 9999/500; * 5% sample;
  ct=1;

* data cleaning;
  nonfilerHH=0; if (txunits=0) then nonfilerHH=1;
  roll = rollovers; if (roll>gross_ret) then roll=gross_ret; * top-code rollovers at gross_ret;
  if (finf_misc > 10000) and (cpftloss=0) then finf_misc = 0; * if a filer had a 1099-MISC income over $10K and no sch C income them set misc income to zero;
  filer=0;  if (tot_inc ne 0) then filer = 1;
  if (misc >= 100000)                then misc = 0; * if a non-filer had a 1099-MISC income of $100K or more then set to zero as they should have filed, so probably problem of uncleaned data;
  miscadj = 0; if (misc>0) then miscadj = 0.3*misc; * Reduce nonfiler 1099-MISC income by 30% because did not deduct losses (only see deductions for filers);

* Income definitions to conform to Census income definition: Personal income plus Social Security, disability and unemployment benefits and less retirement income;
  inc3 = tot_inc + taxexint - dpftloss - finf_ui - txsocsec - txret + wages + div + intr + misc;  * Private non-labor income= 1040 total income + tax exempt interest – Sch D cap gains – filer UI - taxable Soc Sec + nf wages + nf dividends + nf interest + nf misc;  
  inc5 = inc3 + ssa + ui + finf_ui; * Private income plus SS/DI/UI cash transfers = Private income + SS/DI cash transfers + nonfiler UI + filer UI ;
  inc5adj = inc5 - misc + miscadj; * Adjust for misssing non-filer self-employment deductions;

* Income definitions for FULL tax-based income and after-tax income;
  incnew= tot_inc + taxexint - txsocsec + ssa - txret + gross_ret - roll + wages + div + intr + miscadj + ui - dpftloss; * Pre-tax Income;
  inctx = tot_inc + taxexint - txsocsec + ssa - txret + gross_ret - roll + wages + div + intr + miscadj + ui ;  
  incat = inctx + eic + addlchc - totximf; * Federal after-tax Private income plus SS/DI/UI = Private income plus SS/DI/UI cash transfers + Refundable EIC/CTC – Federal Income Tax (no payroll taxes, non-filers have zero income tax);
  incat_norc = inctx - totximf; * After-tax but no ref. credits;
  incnew_noret = incnew - (gross_ret - roll); * Remove private retirement income for CPS comparison;

* Filer income from tax returns and SSA-1099 (no retirement income because do not have info return breakout for filers and generally remove in first half of the paper);
  txret_inc = .; if (tot_inc ne 0) then txret_inc = tot_inc + taxexint - dpftloss - txsocsec - txret + finf_ssa; 

* Size-adjust income definitions and bottom-code at $1;
  szadj_new   = incnew;    if (HHsize>0) then szadj_new   = incnew/sqrt(HHsize);    if (szadj_new  <1)and(szadj_new   ne .) then szadj_new  =1+ranuni(364); 
  szadj_i5    = inc5;      if (HHsize>0) then szadj_i5    = inc5/sqrt(HHsize);      if (szadj_i5   <1)and(szadj_i5    ne .) then szadj_i5   =1+ranuni(364); 
  szadj_i5adj = inc5adj;   if (HHsize>0) then szadj_i5adj = inc5adj/sqrt(HHsize);   if (szadj_i5adj<1)and(szadj_i5adj ne .) then szadj_i5adj=1+ranuni(364); 
  szadj_tx    = inctx;     if (HHsize>0) then szadj_tx    = inctx/sqrt(HHsize);     if (szadj_tx   <1)and(szadj_tx    ne .) then szadj_tx   =1+ranuni(364); 
  szadj_at    = incat;     if (HHsize>0) then szadj_at    = incat/sqrt(HHsize);     if (szadj_at   <1)and(szadj_at    ne .) then szadj_at   =1+ranuni(364); 
  szadj_txr   = txret_inc; if (HHsize>0) then szadj_txr   = txret_inc/sqrt(HHsize); if (szadj_txr  <1)and(szadj_txr   ne .) then szadj_txr  =1+ranuni(364); 
  szadj_nort=incnew_noret; if (HHsize>0) then szadj_nort  = incnew_noret/sqrt(HHsize); if (szadj_nort  <1)and(szadj_nort   ne .) then szadj_nort  =1+ranuni(364); 

  if (inc5adj  <1)and(inc5adj   ne .)   then inc5adj    =1+ranuni(364);
  if (incnew   <1)and(incnew    ne .)   then incnew     =1+ranuni(364);
  if (inctx    <1)and(inctx     ne .)   then inctx      =1+ranuni(364);
  if (incat    <1)and(incat     ne .)   then incat      =1+ranuni(364);
  if (incat_norc<1)and(incat_norc ne .) then incat_norc =1+ranuni(364);
  if (txret_inc<1)and(txret_inc ne .)   then txret_inc  =1+ranuni(364);
  if (incnew_noret<1)and(incnew_noret ne .) then incnew_noret=1+ranuni(364);

  rc=.;  if (eic + addlchc)>0 then rc =eic + addlchc;
  rc2=0; if (eic + addlchc)>0 then rc2=eic + addlchc;
  refcr = 0; if (rc ne .) then refcr=1;
  incat_rc =.; if (refcr=1) then incat_rc =incat;
  incat_rcz=.; if (refcr=1) then incat_rcz=incat_norc;

  eic_cr = 0; if (eic ne .)and(eic > 0) then eic_cr=1;    
  inctaxamt = inctx - incat; 

  txu1=0; if (txunits=1) then txu1=1;
  txu2=0; if (txunits=2) then txu2=1;
  txu3=0; if (txunits=3) then txu3=1;
  txu4=0; if (txunits>3) then txu4=1;

  * Table 3: Number Filing Tax units and Nonfiling Individuals in each household;
  f0_nf1=0; if (txunits-txunits_nf)=0   and txunits_nf=1   then f0_nf1=1; 
  f0_nf2=0; if (txunits-txunits_nf)=0   and txunits_nf>=2  then f0_nf2=1; 
  f1_nf0=0; if (txunits-txunits_nf)=1   and txunits_nf=0   then f1_nf0=1;
  f1_nf1=0; if (txunits-txunits_nf)=1   and txunits_nf=1   then f1_nf1=1;
  f1_nf2=0; if (txunits-txunits_nf)=1   and txunits_nf>=2  then f1_nf2=1; 
  f2_nf1=0; if (txunits-txunits_nf)>=2  and txunits_nf=0   then f2_nf0=1;
  f2_nf2=0; if (txunits-txunits_nf)>=2  and txunits_nf=1   then f2_nf1=1;
  f2_nf3=0; if (txunits-txunits_nf)>=2  and txunits_nf>=2  then f2_nf2=1;

  * Filer income from information returns;
  finf_inc = .;      if (tot_inc ne 0) then finf_inc       = finf_wages + finf_un + finf_ui + finf_misc     + finf_ssa; * W2 wages, 1099-DIV dividends, 1099-INT interest, 1099-G unemployment benefits, SSA-1099 benefits, 1099-MISC;
  adj_finf_inc = .;  if (tot_inc ne 0) then adj_finf_inc   = finf_wages + finf_un + finf_ui + finf_misc*0.3 + finf_ssa; * adjust for avg. net SE income;
  adj_finf_inc50= .; if (tot_inc ne 0) then adj_finf_inc50 = finf_wages + finf_un + finf_ui + finf_misc*0.5 + finf_ssa; * adjust for avg. net SE income;
  finf_misc100  =finf_misc; if (finf_misc100>100000) then finf_misc100=100000;
  finf_misc10050=finf_misc; if (finf_misc10050>(0.5*cpftloss))and(finf_misc100>100000) then finf_misc10050=100000;
  fmisc_cp=.;    if (finf_misc    > 1000) then fmisc_cp   =cpftloss/finf_misc;
  fmisc100_cp=.; if (finf_misc100 > 1000) then fmisc100_cp=cpftloss/finf_misc100;
  if (adj_finf_inc<1)and(adj_finf_inc ne .) then adj_finf_inc = 1+ranuni(364); 
run; 

* Remove group quarters for hhs2, in 2010 this removes 8-9M individuals and about 0.5M households;
data hhs2; set hhs;
   if (HHsize>=11) then delete; 
   ret = 0; taxwt=0; 
   HHwt=HHsize*wt;
run;

* Footnote (NO RETIREMENT INCOME);
%centiles(hhs2, incnew_noret, wt, cent30nr);
proc sort data=hhs2; by cent30nr; run; 
proc means noprint data=hhs2; by cent30nr; var incnew_noret; weight wt;
         output out=Centiles_No_Retire mean(incnew_noret)=mn_incnew; run; 

* Table 2: Number of filing and non-filing tax-units in each household, 2010;
proc means noprint data=hhs (where=(HHsize<11)); var txunits; weight wt;
         output out=TxUnitHHsizes  sum(ct HHsize txunits txunits_nf nonfilerHH f0_nf1 f0_nf2 f1_nf0 f1_nf1 f1_nf2 f2_nf0 f2_nf1 f2_nf2)=
 									   ct HHsize txunits txunits_nf nonfilerHH f0_nf1 f0_nf2 f1_nf0 f1_nf1 f1_nf2 f2_nf0 f2_nf1 f2_nf2; run;

* Table 1: Number of households by HH size, 2010;
proc sort data=hhs; by HHsize; run;
proc means noprint data=hhs (where=(HHsize<11)); by HHsize; var txunits; weight wt; 
         output out=TableHHsize11 sum(ct)=ct; run;

* Not dropping group quarters (HH size of 11 or more); 
proc means noprint data=hhs; by HHsize; var txunits; weight wt; 
         output out=TableHHsizeAll sum(ct)=ct; run;


* Table: Number of individuals and households by state, 2010;
proc sort data=hhs; by state; run;

proc means noprint data=hhs; by state; var txunits; weight wt; 
         output out=IndivsbyState  sum(HHsize)=Hhsize; run;

proc means noprint data=hhs (where=(HHsize<11)); by state; var txunits; weight wt; 
         output out=HHsbyState  sum(HHsize gte20yr txunits_nf ct)=NumIndivs NumAdults NumTaxUnits NumHHs; run;


* Table Appendix: Total income by source, 2010;
proc means noprint data=hhs; var txunits; weight wt;
         output out=HHsumsAll  sum(ct inc5 inc5adj inctx incat HHsize txunits txunits_nf nonfilerHH tot_inc dpftloss taxexint txsocsec txret wage_amt fpftloss cpftloss epftloss ui ssa rollovers gross_ret wages div intr misc finf_ui      totximf eic addlchc )=
 								   ct inc5 inc5adj inctx incat HHsize txunits txunits_nf nonfilerHH tot_inc dpftloss taxexint txsocsec txret wage_amt fpftloss cpftloss epftloss ui ssa rollovers gross_ret wages div intr misc unempcmpinAGI totximf eic addlchc ; run;


* Figure 1: compare filer incomes from tax returns or information returns (order by household tax return based income);
%centiles(hhs2, txret_inc, wt, cent0); * Group into income centiles grouped by number of households;
proc sort data=hhs2; by cent0; run; 
proc means noprint data=hhs2; by cent0; var inc5adj; weight wt;
         output out=filers_infret_byTxRet       
         sum(ct)=TaxUnits mean(txret_inc finf_inc adj_finf_inc adj_finf_inc50 fmisc_cp fmisc100_cp)=
                           mntxret_inc mnfinf_inc adj_finf_inc adj_finf_inc50 fmisc_cp fmisc100_cp
                         sum(cpftloss finf_misc finf_misc100 finf_misc10050)=cpftloss finf_misc finf_misc100 finf_misc10050; run;


* Size-Adjusted Income grouped by #individuals (HHwt);
%centiles(hhs2, szadj_new, HHwt, cent30);
%grp995(hhs2,   szadj_new, HHwt, p995, p995_cutoff);
%grp999(hhs2,   szadj_new, HHwt, p999, p999_cutoff);

* Figures 5,6,7: SIZE-ADJUSTED INCOME;
proc sort data=hhs2; by cent30; run; 
proc means noprint data=hhs2; by cent30; var szadj_new; weight HHwt;
         output out=LorenzNoSOI_SzAdj sum(ct)=HHs mean(szadj_new)=mn_incnew mean(szadj_new)=mn_szadj sum(szadj_new)=sum_incnew
                                      mean(eic_cr)=eic_cr; run; 
* Figure 8: SIZE-ADJUSTED INCOME;
proc sort data=hhs2; by cent30; run; 
proc means noprint data=hhs2; by cent30; var szadj_new; weight wt;
         output out=LorenzNoSOI_EIC sum(ct)=HHs mean(eic)=eic; run; 

* Figure 9: Average tax rates by pre-tax income (includes cap gains);
%centiles(hhs2, szadj_tx, HHwt, cent45);
proc sort data=hhs2; by cent45; run; 
proc means noprint data=hhs2; by cent45; var incnew; weight wt;
         output out=TaxRates sum(ct)=HHs mean(inctx incat inctaxamt)=inctx incat inctaxamt; run; 
 
* Table 3 (size-adjust): Pre-tax Size-Adjusted Income percentiles and groups (Census type income + private retirement income);
data hhs2; set hhs2; 
q1=.; q2=.; q3=.; q4=.; top8095=.; top5=.;
       if (cent30<=20) then q1=szadj_new; *incnew; 
  else if (cent30<=40) then q2=szadj_new; *incnew; 
  else if (cent30<=60) then q3=szadj_new; *incnew; 
  else if (cent30<=80) then q4=szadj_new; *incnew; 
  else if (cent30<=95) then top8095=szadj_new; *incnew; 
  else if (cent30 >95) then top5=szadj_new; *incnew;
  pos5=szadj_new; if (szadj_new<0) then pos5=0;
  q1pos=q1; if (szadj_new<0) then q1pos=0;
  top1=.; if (cent30=100) then top1=szadj_new; *incnew;
  year=&yr;
run;

proc means noprint data=hhs2; var incnew; weight HHwt; * use indiv wts;
         output out=Ineq_SzAdj mean(year)=year sum(ct szadj_new pos5 q1 q1pos q2 q3 q4 top8095 top5)=ct tot_incnew tot_pos q1 q1pos q2 q3 q4 top8095 top5
                             P10(szadj_new)=p10 p20(szadj_new)=p20 p40(szadj_new)=p40 p50(szadj_new)=p50 p60(szadj_new)=p60 p80(szadj_new)=p80 p90(szadj_new)=p90 p95(szadj_new)=p95 p99(szadj_new)=p99
                             mean(q1 q1pos q2 q3 q4  top8095 top5)=mn_q1 mn_q1pos mn_q2 mn_q3 mn_q4  mn_top8095 mn_top5 
                             sum(top1 p995 p999)=top1 top05 top01 min(p995 p999)=p995 p999 mean(top1 p995 p999)=mn_top1 mn_top05 mn_top01; run;

%gini3(szadj_new, HHwt, hhs2, &yr, _new);


* Not size-adjusted;
%centiles(hhs2, incnew, wt, cent30);
%grp995(hhs2,   incnew, wt, p995, p995_cutoff);
%grp999(hhs2,   incnew, wt, p999, p999_cutoff);

proc sort data=hhs2; by cent30; run; 
proc means noprint data=hhs2; by cent30; var incnew; weight wt;
         output out=LorenzNoSOI_New sum(ct)=HHs mean(incnew)=mn_incnew sum(incnew)=sum_incnew
         mean(eic_cr eic txunits txu1 txu2 txu3 txu4)=eic_cr eic mn_txunits txu1 txu2 txu3 txu4 ; run; 

* Figure 9: Average tax rates by pre-tax income (includes cap gains);
%centiles(hhs2, inctx, wt, cent40);
proc sort data=hhs2; by cent40; run; 
proc means noprint data=hhs2; by cent40; var incnew; weight wt;
         output out=TaxRates sum(ct)=HHs mean(inctx incat inctaxamt)=inctx incat inctaxamt; run; 

* Table 4 (No SOI): Pre-tax Income percentiles and groups (Census type income + private retirement income);
data hhs2; set hhs2;
q1=.; q2=.; q3=.; q4=.; top8095=.; top5=.;
       if (cent30<=20) then q1=incnew; 
  else if (cent30<=40) then q2=incnew; 
  else if (cent30<=60) then q3=incnew; 
  else if (cent30<=80) then q4=incnew; 
  else if (cent30<=95) then top8095=incnew; 
  else if (cent30 >95) then top5=incnew;
  pos5=incnew; if (incnew<0) then pos5=0;
  q1pos=q1; if (incnew<0) then q1pos=0;
  top1=.; if (cent30=100) then top1=incnew;
  year=&yr;
run;

proc means noprint data=hhs2; var incnew; weight wt;
         output out=Ineq_NoSOI_New mean(year)=year sum(ct incnew pos5 q1 q1pos q2 q3 q4 top8095 top5)=ct tot_incnew tot_pos q1 q1pos q2 q3 q4 top8095 top5
                             P10(incnew)=p10 p20(incnew)=p20 p40(incnew)=p40 p50(incnew)=p50 p60(incnew)=p60 p80(incnew)=p80 p90(incnew)=p90 p95(incnew)=p95 p99(incnew)=p99
                             mean(q1 q1pos q2 q3 q4  top8095 top5)=mn_q1 mn_q1pos mn_q2 mn_q3 mn_q4  mn_top8095 mn_top5 
                             sum(top1 p995 p999)=top1 top05 top01 min(p995 p999)=p995 p999 mean(top1 p995 p999)=mn_top1 mn_top05 mn_top01; run;

%gini3(incnew, wt, hhs2, &yr, _new);


%mend CDWdata;


libname ITM '...';
libname INDIV '...';

%let year = 2010;


%CDWdata(&year);
%TaxUnitDist2010;