%ACCOUNTS    Accounting used for 
%
%              Technology Capital and the Current Account
%              by Ellen McGrattan and Edward Prescott
%

%            Ellen McGrattan, 2-28-08


%
% Sample t0:t1
%

t0           = 1960;
t1           = 2006;

%
% NIPA Table 1.1.5 Gross domestic product (Annual, Bil$)
%

load bea/nipa/table1_1_5A.dat
[n,m]        = size(table1_1_5A);
t            = table1_1_5A(1,:);
i            = find(t==t0):find(t==t1);
table1_1_5A  = table1_1_5A(2:n,i);

gdp          = table1_1_5A(1,:)';
pce          = table1_1_5A(2,:)';
pced         = table1_1_5A(3,:)';
pcend        = table1_1_5A(4,:)';
pces         = table1_1_5A(5,:)';
gpdi         = table1_1_5A(6,:)';
gpdif        = table1_1_5A(7,:)';
gpdii        = table1_1_5A(12,:)';
netxm        = table1_1_5A(13,:)';
netx         = table1_1_5A(14,:)';
netm         = table1_1_5A(17,:)';
gov          = table1_1_5A(20,:)';

%
% NIPA Table 1.1.9 Implicit price deflators (2000=100)
%

load bea/nipa/table1_1_9A.dat
[n,m]        = size(table1_1_9A);
t            = table1_1_9A(1,:);
i            = find(t==t0):find(t==t1);
table1_1_9A  = table1_1_9A(2:n,i);

dgdp0        = table1_1_9A(1,:)';
dgdp         = dgdp0/dgdp0(1);

%
% NIPA Table 1.7.5 Relation of GDP, GNP, etc. (Annual, Bil$)
%

load bea/nipa/table1_7_5A.dat
[n,m]        = size(table1_7_5A);
t            = table1_7_5A(1,:);
i            = find(t==t0):find(t==t1);
table1_7_5A  = table1_7_5A(2:n,i);

gnp          = table1_7_5A(4,:)';
gni          = table1_7_5A(28,:)';

%
% NIPA Table 1.10 Gross domestic income (Annual, Bil$)
%

load bea/nipa/table1_10A.dat
[n,m]        = size(table1_10A);
t            = table1_10A(1,:);
i            = find(t==t0):find(t==t1);
table1_10A   = table1_10A(2:n,i);

gdi          = table1_10A(1,:)';
comp         = table1_10A(2,:)';
ibt          = table1_10A(9,:)'-table1_10A(10,:)'+table1_10A(14,:)';
nint         = table1_10A(13,:)';
pinc         = table1_10A(15,:)';
rent         = table1_10A(16,:)';
prof         = table1_10A(17,:)';
gprof        = table1_10A(22,:)';
dep          = table1_10A(23,:)';
desc         = table1_10A(26,:)';

%
% NIPA Table 1.12  National income by type of income (Anuual, Bil$)
%

load bea/nipa/table1_12A.dat
[n,m]        = size(table1_12A);
t            = table1_12A(1,:);
i            = find(t==t0):find(t==t1);
table1_12A   = table1_12A(2:n,i);

natprof      = table1_12A(13,:)';

%
% NIPA Table 1.13  National income by sector (Anuual, Bil$)
%

load bea/nipa/table1_13A.dat
[n,m]        = size(table1_13A);
t            = table1_13A(1,:);
i            = find(t==t0):find(t==t1);
table1_13A   = table1_13A(2:n,i);

dbus         = table1_13A(2,:)';
dcbus        = table1_13A(3,:)';
dcbuscomp    = table1_13A(4,:)';
dcbusprof    = table1_13A(7,:)';
dcbusnint    = table1_13A(8,:)';
dcbusibt     = table1_13A(9,:)';
dncbus       = table1_13A(10,:)';
dncbuscomp   = table1_13A(11,:)';
dncbuspinc   = table1_13A(14,:)';
dncbusrent   = table1_13A(15,:)';
dncbusnint   = table1_13A(16,:)';
dncbusibt    = table1_13A(17,:)';
dncbusgprof  = table1_13A(18,:)';
dppbus       = table1_13A(19,:)';
dppbuscomp   = table1_13A(20,:)';
dppbuspinc   = table1_13A(23,:)';
dppbusfpinc  = table1_13A(24,:)';
dppbusnfpinc = table1_13A(25,:)';
dppbusnint   = table1_13A(26,:)';
dppbusibt    = table1_13A(27,:)';
dopbus       = table1_13A(28,:)';
dopbuscomp   = table1_13A(29,:)';
dopbuspinc   = table1_13A(32,:)';
dopbusrent   = table1_13A(33,:)';
dopbusnint   = table1_13A(34,:)';
dopbusibt    = table1_13A(35,:)';
dgebus       = table1_13A(36,:)';
dgebuscomp   = table1_13A(37,:)';
dgebusgprof  = table1_13A(40,:)';
dhh          = table1_13A(42,:)';
dhhcomp      = table1_13A(43,:)';
dhhrent      = table1_13A(46,:)';
dhhnint      = table1_13A(47,:)';
dhhibt       = table1_13A(48,:)';
dnp          = table1_13A(49,:)';
dnpcomp      = table1_13A(50,:)';
dnprent      = table1_13A(53,:)';
dnpnint      = table1_13A(54,:)';
dnpibt       = table1_13A(55,:)';
dgg          = table1_13A(56,:)';
dggcomp      = table1_13A(57,:)';
row          = table1_13A(60,:)';
rowcomp      = table1_13A(61,:)';
rowprof      = table1_13A(62,:)';
rownint      = table1_13A(63,:)';

ssectcomp    = dcbuscomp+dppbuscomp+dopbuscomp;
ssectprof    = dcbusprof;
ssectrent    = dopbusrent;
ssectpinc    = dppbuspinc+dopbuspinc;
ssectnint    = dcbusnint+dppbusnint+dopbusnint;
ssectibt     = dcbusibt+dppbusibt+dopbusibt;
ssecti       = ssectcomp+ssectprof+ssectrent+ssectpinc+ssectnint+ssectibt;

osectcomp    = dhhcomp+dnpcomp+dgebuscomp+dggcomp;
osectrent    = dhhrent+dnprent;
osectprof    = dgebusgprof;
osectnint    = dhhnint+dnpnint;
osectibt     = dhhibt+dnpibt;
osecti       = osectcomp+osectprof+osectrent+osectnint+osectibt;

%
% NIPA Table 1.14  Gross VA of domestic corporate business (Annual, Bil$)
%

load bea/nipa/table1_14A.dat
[n,m]        = size(table1_14A);
t            = table1_14A(1,:);
i            = find(t==t0):find(t==t1);
table1_14A   = table1_14A(2:n,i);

cgdi         = table1_14A(1,:)';
cdep         = table1_14A(2,:)';
ccomp        = table1_14A(4,:)';
cibt         = table1_14A(7,:)'+table1_14A(10,:)';
cnint        = table1_14A(9,:)';
cprof        = table1_14A(11,:)';
ctax         = table1_14A(12,:)';
ndiv         = table1_14A(14,:)';
rete         = table1_14A(15,:)';
pbtd         = table1_14A(32,:)';
patd         = table1_14A(33,:)';
ivadj        = table1_14A(34,:)';  
ccadj        = table1_14A(35,:)';  

%
% NIPA Table 2.5.5 Personal consumption expenditures by type (Annual, Bil$)
%

load bea/nipa/table2_5_5A.dat
[n,m]        = size(table2_5_5A);
t            = table2_5_5A(1,:);
i            = find(t==t0):find(t==t1);
table2_5_5A  = table2_5_5A(2:n,i);

pbexp        = table2_5_5A(64,:)';

%
% NIPA Table 3.1 Govt current receipts and expenditures (Annual, Bil$)
%

load bea/nipa/table3_1A.dat
[n,m]        = size(table3_1A);
t            = table3_1A(1,:);
i            = find(t==t0):find(t==t1);
table3_1A    = table3_1A(2:n,i);

govc         = table3_1A(16,:)';
govi         = table3_1A(35,:)';
hhdt         = table3_1A(3,:)';
hhss         = table3_1A(7,:)';


%
% NIPA Table 3.2 Federal govt current receipts and expenditures (Annual, Bil$)
%

load bea/nipa/table3_2A.dat
[n,m]        = size(table3_2A);
t            = table3_2A(1,:);
i            = find(t==t0):find(t==t1);
table3_2A    = table3_2A(2:n,i);

fedres       = table3_2A(8,:)';

%
% NIPA Table 3.5 Taxes on production and imports (Annual, Bil$)
%

load bea/nipa/table3_5A.dat
[n,m]        = size(table3_5A);
t            = table3_5A(1,:);
i            = find(t==t0):find(t==t1);
table3_5A    = table3_5A(2:n,i);

ibtonly      = table3_5A(1,:)';
stax         = sum(table3_5A([3,11,14,28,29,30,31],:))';
ptax         = table3_5A(27,:)';

%
% NIPA Table 4.1A Foreign Transactions in NIPA (Annual, Bil$)
%

load bea/nipa/table4_1A.dat
[n,m]        = size(table4_1A);
t            = table4_1A(1,:);
i            = find(t==t0):find(t==t1);
table4_1A    = table4_1A(2:n,i);

ex           = table4_1A(2,:)';
im           = table4_1A(14,:)';
nx           = ex-im;
nfrinc       = table4_1A(7,:)';
nfrcomp      = table4_1A(8,:)';
nfrint       = table4_1A(10,:)';
nfrdiv       = table4_1A(11,:)';
nfrre        = table4_1A(12,:)';
nfrdipi      = nfrint+nfrdiv+nfrre;
nfpinc       = table4_1A(19,:)';
nfpcomp      = table4_1A(20,:)';
nfpint       = table4_1A(22,:)';
nfpdiv       = table4_1A(23,:)';
nfpre        = table4_1A(24,:)';
nfpdipi      = nfpint+nfpdiv+nfpre;
tran         = table4_1A(25,:)';
ca           = table4_1A(29,:)';

%
% NIPA Table 5.7.5A Private inventories and final sales (Quarterly, Bil$)
%

load bea/nipa/table5_7_5BQ.dat
load bea/nipa/table5_7_5AQ.dat
table5_7_5BQ = [table5_7_5BQ(1,:)+table5_7_5BQ(2,:)/10;
                table5_7_5BQ(3:4,:)];
table5_7_5AQ = [table5_7_5AQ(1,:)+table5_7_5AQ(2,:)/10;
                table5_7_5AQ(3:4,:)];
i            = find(table5_7_5AQ(1,:)==table5_7_5BQ(1,1));
table5_7_5Q  = [table5_7_5AQ(1:3,1:i-1),table5_7_5BQ(1:3,:)];
table5_7_5A  = (qtr2ann(table5_7_5Q',1))';
[n,m]        = size(table5_7_5A);
t            = table5_7_5A(1,:);
i            = find(t==t0):find(t==t1);
table5_7_5A  = table5_7_5A(2:n,i);

kcinv        = table5_7_5A(1,:)'-.9*table5_7_5A(2,:)';
kncinv       = .9*table5_7_5A(2,:)';

%
% NIPA Table 7.5 Consumption of fixed capital by legal form (Annual, Bil$)
%

load bea/nipa/table7_5A.dat
[n,m]        = size(table7_5A);
t            = table7_5A(1,:);
i            = find(t==t0):find(t==t1);
table7_5A    = table7_5A(2:n,i);

dbusdep      = table7_5A(3,:)';
dcbusdep     = table7_5A(4,:)';
dncbusdep    = table7_5A(7,:)';
dppbusdep    = table7_5A(8,:)';
dopbusdep    = table7_5A(11,:)';
dgebusdep    = table7_5A(25,:)';
dhhdep       = table7_5A(19,:)';
dnpdep       = table7_5A(20,:)';
dggdep       = table7_5A(22,:)';

ssectdep     = dcbusdep+dppbusdep+dopbusdep;
ssect        = ssecti+ssectdep;
osectdep     = dhhdep+dnpdep+dgebusdep+dggdep;
osect        = osecti+osectdep;

%
% NIPA Table 7.19  Relation of corporate profits, etc. in NIPA and IRS
%

load bea/nipa/table7_19A.dat
[n,m]        = size(table7_19A);
table7_19A   = [table7_19A,[2002:2004;NaN*ones(n-1,3)]];
t            = table7_19A(1,:);
i            = find(t==t0):find(t==t1);
table7_19A   = table7_19A(2:n,i);

gain         = table7_19A(12,:)';

%
% FA Table 1.1 Current-cost net stocks (Annual, Bil$, year-end) 
%

load bea/fa/table1_1A.dat
[n,m]        = size(table1_1A);
t            = table1_1A(1,:);
i            = find(t==t0):find(t==t1);
table1_1A    = table1_1A(2:n,i);

kfix         = table1_1A(2,:)';
kgfix        = table1_1A(8,:)';
kcd          = table1_1A(13,:)';

%
% FA Table 1.5 Historical-cost investments (Annual, Bil$, year-end) 
%

load bea/fa/table1_5A.dat
[n,m]        = size(table1_5A);
t            = table1_5A(1,:);
i            = find(t==t0):find(t==t1);
table1_5A    = table1_5A(2:n,i);

xfix         = table1_5A(2,:)';
xgfix        = table1_5A(8,:)';
xcd          = table1_5A(13,:)';

%
% FA Tables 6.1 Current-cost net stocks by owner (Annual, Bil$, year-end) 
%

load bea/fa/table6_1A.dat
[n,m]        = size(table6_1A);
t            = table6_1A(1,:);
i            = find(t==t0):find(t==t1);
t6_1A        = table6_1A(1,i)';
table6_1A    = table6_1A(2:n,i);

kcfix        = table6_1A(2,:)';

kcnffix        = table6_1A(4,:)';

kppfix       = table6_1A(6,:)'+table6_1A(7,:)';
khhfix       = table6_1A(9,:)';
knpfix       = table6_1A(8,:)';
ktefix       = table6_1A(10,:)';
kncfix       = kppfix+khhfix+knpfix+ktefix;

%
% FA Table 6.7 Historical-cost investments by owner (Annual, Bil$, year-end) 
%

load bea/fa/table6_7A.dat
[n,m]        = size(table6_7A);
t            = table6_7A(1,:);
i            = find(t==t0):find(t==t1);
t6_7A        = table6_7A(1,i)';
table6_7A    = table6_7A(2:n,i);

xcfix        = table6_7A(2,:)';
xppfix       = table6_7A(6,:)'+table6_7A(7,:)';
xhhfix       = table6_7A(9,:)';
xnpfix       = table6_7A(8,:)';
xtefix       = table6_7A(10,:)';
xncfix       = xppfix+xhhfix+xnpfix+xtefix;

%
% FOF Table F6 Gross domestic product
%

load fof/atabs/atab6d.dat
atab6dA      = (qtr2ann(atab6d))';
[n,m]        = size(atab6dA);
t            = atab6dA(1,:);
i            = find(t==t0):find(t==t1);
atab6dA      = atab6dA(2:n,i)/1000;

cgpdif       = sum(atab6dA([10,13,16,18],:))'+.1*atab6dA(12,:)';
npgpdif      = atab6dA(9,:)';
hhgpdif      = atab6dA(15,:)';
ncgpdif      = sum(atab6dA([11,17],:))'+.9*atab6dA(12,:)';

cgpdii       = atab6dA(20,:)'+.1*atab6dA(22,:)';
ncgpdii      = atab6dA(21,:)'+.9*atab6dA(22,:)';


%
% FOF Table F10 Measures of personal saving (Quarterly, Mil$)
%
load fof/atabs/atab10d.dat
atab10dA     = (qtr2ann(atab10d))';
[n,m]        = size(atab10dA);
t            = atab10dA(1,:);
i            = find(t==t0):find(t==t1);
atab10dA     = atab10dA(2:n,i)/1000;

cddep        = atab10dA(27,:)';

%
% FOF Table B100 Balance sheet of households (Quarterly, Bil$)
%

load fof/btabs/btab100d.dat
btab100dA    = (qtr2ann(btab100d,1))';
[n,m]        = size(btab100dA);
t            = btab100dA(1,:);
i            = find(t==t0):find(t==t1);
btab100dA    = btab100dA(2:n,i);

hhland       = btab100dA(3,:)'-btab100dA(42,:)'-btab100dA(46,:)';
eqnc         = btab100dA(28,:)';

%
% FOF Table L213 Corporate Equities (Quarterly, Bil$)
%

load fof/ltabs/ltab213d.dat
ltab213dA    = (qtr2ann(ltab213d,1))';
[n,m]        = size(ltab213dA);
t            = ltab213dA(1,:);
i            = find(t==t0):find(t==t1);
ltab213dA    = ltab213dA(2:n,i)/1000;

eqc          = ltab213dA(20,:)';

%
% FOF Table B102 Balance sheet of nonfarm nonfin. corporate (Quarterly, Bil$)
%

load fof/btabs/btab102d.dat
btab102dA    = (qtr2ann(btab102d,1))';
[n,m]        = size(btab102dA);
t            = btab102dA(1,:);
i            = find(t==t0):find(t==t1);
btab102dA    = btab102dA(2:n,i);

cland        = btab102dA(3,:)'-btab102dA(33,:)'-btab102dA(34,:)';

%
% FOF Table B103 Balance sheet of nonfarm noncorporate (Quarterly, Bil$)
%

load fof/btabs/btab103d.dat
btab103dA    = (qtr2ann(btab103d,1))';
[n,m]        = size(btab103dA);
t            = btab103dA(1,:);
i            = find(t==t0):find(t==t1);
btab103dA    = btab103dA(2:n,i);

nfbland      = btab103dA(3,:)'-btab103dA(32,:)'-btab103dA(33,:)';

land         = hhland+cland+nfbland;
k            = kfix+kcinv+kncinv+kcd+land;

%
% FOF Table L107 Rest of World (Quarterly, Bil$)
%

load fof/ltabs/ltab107d.dat
ltab107dA    = (qtr2ann(ltab107d,1))';
[n,m]        = size(ltab107dA);
t            = ltab107dA(1,:);
i            = find(t==t0):find(t==t1);
ltab107dA    = ltab107dA(2:n,i)/1000;

usdia_cc     = ltab107dA(42,:)';
fdius_cc     = ltab107dA(23,:)';


%
% ITA Table 1 US International Transactions (Annual, Mil$)
%
load bea/ita/Table1A.dat
[n,m]        = size(Table1A);
t            = Table1A(1,:);
i            = find(t==t0):find(t==t1);
ita1         = Table1A(2:n,i);

usdia        = ita1(14,:)'/1000;
fdius        = -ita1(31,:)'/1000;

%
% ITA Table 6 Direct Investment Incomes (Annual, Mil$)
%
load bea/ita/Table6aA.dat
load bea/ita/Table6bA.dat
load bea/scb/DIAint.dat
[n,m]        = size(Table6aA);
t            = 1957:Table6aA(1,m);
usdia_int    = [DIAint(:,2)',Table6bA(6,1:17),Table6aA(8,:)];
fdius_int    = -[DIAint(:,3)',Table6bA(47,1:17),Table6aA(73,:)];
i            = find(t==t0):find(t==t1);
usdia_int    = usdia_int(:,i)'/1000;
fdius_int    = fdius_int(:,i)'/1000;
dianint      = usdia_int-fdius_int;

%
% IIP Table 2 Direct Investment Positions (Annual, Mil$)
%

load bea/iip/dipos.dat
[n,m]        = size(dipos);
t            = dipos(:,1);
i            = find(t==t0):find(t==t1);
dipos        = dipos(i,2:m)/1000;

usdia_cc2    = dipos(:,1);
fdius_cc2    = dipos(:,3);
usdia_mv     = dipos(:,2);
fdius_mv     = dipos(:,4);


%
% Taxes
%
i            = .041;
mu1          = (cibt+dppbusibt+dopbusibt)./ibt;
mu2          = 1-mu1;
tem          = stax./(pce+i*kcd+cddep);
tauc         = tem./(1-tem);
cdtax        = tauc.*pced./(1+tauc);
cntax        = stax-cdtax; 
taup         = (ctax-fedres)./(cprof-fedres);
tauh1        = 1.6*hhdt./(gdp-dep-ibtonly)+hhss./(comp+.7*pinc);
tauh2        = 1.6*hhdt./(gain+gdp-dep-ibtonly)+hhss./(comp+.7*pinc);

disp('Normalize variables with:')
disp('  1. Adjusted GNP ')
disp('  2. Adjusted GDP ')
disp('  3. NIPA GNP ')
disp('  4. NIPA GDP ')
inorm        = input('Choice: ');
adj          = -stax-pbexp+cddep+i*kcd+i*kgfix;
agnp         = gnp+adj;
agni         = gni+adj;
agdp         = gdp+adj;
if inorm==1;
 y = agnp;
elseif inorm==2;
 y = agdp;
elseif inorm==3;
 y = gnp;
else
 y = gdp;
end;

clc
t    = [t0:t1]';
fprintf('US ADJUSTED NATIONAL ACCOUNTS, %4.0f:%4.0f\n',t0,t1)
fprintf('============================================\n')
fprintf('INCOME\n')
fprintf(' Business                        %6.3f\n', ...
            mean((ssect-mu1.*stax-pbexp)./y))
fprintf('   Capital consumption           %6.3f\n',mean(ssectdep./y))
fprintf('   Business taxes                %6.3f\n',mean(ssectibt./y))
fprintf('   Less: sales tax               %6.3f\n',mean(mu1.*stax./y))
fprintf('   Compensation                  %6.3f\n',mean(ssectcomp./y))
fprintf('   Proprietors'' income           %6.3f\n',mean(ssectpinc./y))
fprintf('   Corporate profits             %6.3f\n',mean(ssectprof./y))
fprintf('   Rental income                 %6.3f\n',mean(ssectrent./y))
fprintf('   Net interest                  %6.3f\n',mean(ssectnint./y))
fprintf('   Less: intermediate service    %6.3f\n',mean(pbexp./y))
fprintf('\n')
fprintf(' Nonbusiness                     %6.3f\n', ...
            mean((osect+cddep-mu2.*stax+i*kgfix+i*kcd)./y))
fprintf('   Capital consumption           %6.3f\n',mean(osectdep./y))
fprintf('   Plus: durable depreciation    %6.3f\n',mean(cddep./y))
fprintf('   Business taxes                %6.3f\n',mean(osectibt./y))
fprintf('   Less: sales tax               %6.3f\n',mean(mu2.*stax./y))
fprintf('   Compensation                  %6.3f\n',mean(osectcomp./y))
fprintf('   Govt. enterprise profits      %6.3f\n',mean(osectprof./y))
fprintf('   Rental income                 %6.3f\n',mean(osectrent./y))
fprintf('   Net interest                  %6.3f\n',mean(osectnint./y))
fprintf('   Plus: imputed services        %6.3f\n', ...
            i*mean((kgfix+kcd)./y))
fprintf(' --------------------------------------\n')
fprintf(' Adjusted gross domestic income  %6.3f\n',mean((gdi+adj)./y))
fprintf('\n')
fprintf(' Net factor income               %6.3f\n',mean(row./y))
fprintf('   Net compensation              %6.3f\n',mean(rowcomp./y))
fprintf('   Net profits+interest          %6.3f\n',mean((rowprof+rownint)./y))
fprintf('     DI receipts                 %6.3f\n',mean(usdia./y))
fprintf('     Less: DI payments           %6.3f\n',mean(fdius./y))
fprintf('     Other receipts              %6.3f\n',mean((nfrdipi-usdia)./y))
fprintf('     Less: other payments        %6.3f\n',mean((nfpdipi-fdius)./y))
fprintf('============================================\n')
fprintf(' Adjusted gross national income  %6.3f\n',mean((gni+adj)./y))
fprintf('\n')
fprintf('+Statistical descrepancy         %6.3f\n',mean(desc./y))
fprintf('============================================\n')
fprintf('=Adjusted gross national product %6.3f\n',mean((gnp+adj)./y))
fprintf('============================================\n')
fprintf('\n')
fprintf('\n')
fprintf('PRODUCT\n')
fprintf(' Consumption                     %6.3f\n', ...
            mean((pce-pced+i*kcd+cddep-pbexp-cntax+govc+i*kgfix)./y))
fprintf('   PCE nondurable+service        %6.3f\n',mean((pce-pced)./y))
fprintf('   Plus: imputed service         %6.3f\n',mean(i*kcd./y))
fprintf('   Plus: durable depreciation    %6.3f\n',mean(cddep./y))
fprintf('   Less: intermediate service    %6.3f\n',mean(pbexp./y))
fprintf('   Less: sales tax               %6.3f\n',mean(cntax./y))
fprintf('   Government consumption        %6.3f\n',mean(govc./y))
fprintf('   Plus: imputed service         %6.3f\n',mean(i*kgfix./y))
fprintf('\n')
fprintf(' Tangible investment             %6.3f\n', ...
            mean((pced+gpdi+govi-cdtax)./y))
fprintf('   Business                      %6.3f\n', ...
              mean((cgpdif+cgpdii+ncgpdif+ncgpdii)./y))
fprintf('     Corporate                   %6.3f\n',mean((cgpdif+cgpdii)./y))
fprintf('     Noncorporate                %6.3f\n',mean((ncgpdif+ncgpdii)./y))
fprintf('   Nonbusiness                   %6.3f\n', ...
              mean((hhgpdif+npgpdif+govi+pced-cdtax)./y))
fprintf('     Household GPDI              %6.3f\n',mean(hhgpdif./y))
fprintf('     Nonprofits GPDI             %6.3f\n',mean(npgpdif./y))
fprintf('     Government investment       %6.3f\n',mean(govi./y))
fprintf('     Consumer durables           %6.3f\n',mean(pced./y))
fprintf('     Less: sales tax             %6.3f\n',mean(cdtax./y))
fprintf(' Net Exports                     %6.3f\n',mean(netxm./y))
fprintf(' --------------------------------------\n')
fprintf(' Adjusted gross domestic product %6.3f\n',mean((gdp+adj)./y))
fprintf('\n')
fprintf('+Net factor income               %6.3f\n',mean(row./y))
fprintf('============================================\n')
fprintf('=Adjusted gross national product %6.3f\n',mean((gnp+adj)./y))
fprintf('============================================\n')
fprintf(' NOTE: adjustments made are\n')
fprintf('\n')
fprintf('  1. Removed sales tax from consumption and IBT\n')
fprintf('  2. Added private and public consumption services\n')
fprintf('  3. Removed intermediate financial services\n')
fprintf('\n')
fprintf('\n')
fprintf('US TAX RATES, %4.0f:%4.0f\n',t0,t1)
fprintf('============================================\n')
fprintf('   Consumption                   %6.3f\n',mean(tauc))
fprintf('   Labor\n')
fprintf('     Including capital gains     %6.3f\n',mean(tauh1))
fprintf('     Excluding capital gains     %6.3f\n',avgnan(tauh2))
fprintf('   Corporate profits             %6.3f\n',mean(taup))
fprintf('============================================\n')
fprintf('\n')
fprintf('\n')


itab4     = input('Print Table 4? (1=yes, 0=no) ');

if itab4;
fid  = fopen('table4.out','w');
fprintf(fid,'US ADJUSTED NATIONAL ACCOUNTS, %4.0f:%4.0f\n',t0,t1);
fprintf(fid,'============================================\n');
fprintf(fid,'INCOME\n');
fprintf(fid,' Business                        %6.3f\n', ...
                mean((ssect-mu1.*stax-pbexp)./y));
fprintf(fid,'   Depreciation                  %6.3f\n',mean(ssectdep./y));
fprintf(fid,'   Capital income                %6.3f\n', ...
                mean((.3*ssectpinc+ssectprof+ssectrent+ssectnint+ ...
                ssectibt-mu1.*stax-pbexp)./y));
fprintf(fid,'   Labor income                  %6.3f\n', ...
                mean((ssectcomp+.7*ssectpinc)./y));
fprintf(fid,' Nonbusiness                     %6.3f\n', ...
                mean((osect+cddep-mu2.*stax+i*kgfix+i*kcd)./y));
fprintf(fid,' --------------------------------------\n');
fprintf(fid,' Adjusted gross domestic income  %6.3f\n',mean((gdi+adj)./y));
fprintf(fid,' Net factor income               %6.3f\n',mean(row./y));
fprintf(fid,'+Statistical descrepancy         %6.3f\n',mean(desc./y));
fprintf(fid,'============================================\n');
fprintf(fid,'=Adjusted gross national product %6.3f\n',mean((gnp+adj)./y));
fprintf(fid,'============================================\n');
fprintf(fid,'\n');
fprintf(fid,'PRODUCT\n');
fprintf(fid,' Consumption                     %6.3f\n', ...
              mean((pce-pced+i*kcd+cddep-pbexp-cntax+govc+i*kgfix)./y));
fprintf(fid,' Tangible investment             %6.3f\n', ...
              mean((pced+gpdi+govi-cdtax)./y));
fprintf(fid,'   Business                      %6.3f\n', ...
                mean((cgpdif+cgpdii+ncgpdif+ncgpdii)./y));
fprintf(fid,'   Nonbusiness                   %6.3f\n', ...
                mean((hhgpdif+npgpdif+govi+pced-cdtax)./y));
fprintf(fid,' Net Exports                     %6.3f\n',mean(netxm./y));
fprintf(fid,' --------------------------------------\n');
fprintf(fid,' Adjusted gross domestic product %6.3f\n',mean((gdp+adj)./y));
fprintf(fid,'+Net factor income               %6.3f\n',mean(row./y));
fprintf(fid,'   Net compensation              %6.3f\n',mean(rowcomp./y));
fprintf(fid,'   Net profits+interest          %6.3f\n', ...
                  mean((rowprof+rownint)./y));
fprintf(fid,'     DI receipts                 %6.3f\n',mean(usdia./y));
fprintf(fid,'     Less: DI payments           %6.3f\n',mean(fdius./y));
fprintf(fid,'     Other receipts              %6.3f\n', ...
                  mean((nfrdipi-usdia)./y));
fprintf(fid,'     Less: other payments        %6.3f\n', ...
                  mean((nfpdipi-fdius)./y));
fprintf(fid,'============================================\n');
fprintf(fid,'=Adjusted gross national product %6.3f\n',mean((gnp+adj)./y));
fprintf(fid,'============================================\n');
fclose(fid);
end;


iextra = input('Print details? (1=yes, 0=no) ');
if iextra ==1; 
fprintf('\n')
fprintf('\n')
fprintf('DETAILS: NIPA AND FIXED ASSET AVERAGES\n',t0,t1)
fprintf('============================================\n')
fprintf(' Gross Domestic Income           %6.3f\n',mean(gdi./y))
fprintf('   Capital consumption           %6.3f\n',mean(dep./y))
fprintf('   IBT+transfers-subsidies       %6.3f\n',mean(ibt./y))
fprintf('   Compensation                  %6.3f\n',mean(comp./y))
fprintf('   Corporate Profits (adj.)      %6.3f\n',mean(prof./y))
fprintf('   Proprietors'' income           %6.3f\n',mean(pinc./y))
fprintf('   Rental income                 %6.3f\n',mean(rent./y))
fprintf('   Govt. enterprises             %6.3f\n',mean(gprof./y))
fprintf('   Net interest                  %6.3f\n',mean(nint./y))
fprintf('\n')
fprintf('\n')
fprintf(' BY SECTOR:\n')
fprintf('   Domestic corporate            %6.3f\n',mean(cgdi./y))
fprintf('     Capital consumption         %6.3f\n',mean(cdep./y))
fprintf('     IBT+transfers-subsidies     %6.3f\n',mean(cibt./y))
fprintf('     Compensation                %6.3f\n',mean(ccomp./y))
fprintf('     Profits (adj)               %6.3f\n',mean(prof./y))
fprintf('     Net interest                %6.3f\n',mean(cnint./y))
fprintf('\n')
fprintf('   Proprietors/Partners          %6.3f\n',mean((dppbus+dppbusdep)./y))
fprintf('     Capital consumption         %6.3f\n',mean(dppbusdep./y))
fprintf('     IBT+transfers-subsidies     %6.3f\n',mean(dppbusibt./y))
fprintf('     Compensation                %6.3f\n',mean(dppbuscomp./y))
fprintf('     Proprietors'' income         %6.3f\n',mean(dppbuspinc./y))
fprintf('     Net interest                %6.3f\n',mean(dppbusnint./y))
fprintf('\n')
fprintf('   Other Private                 %6.3f\n',mean((dopbus+dopbusdep)./y))
fprintf('     Capital consumption         %6.3f\n',mean(dopbusdep./y))
fprintf('     IBT+transfers-subsidies     %6.3f\n',mean(dopbusibt./y))
fprintf('     Compensation                %6.3f\n',mean(dopbuscomp./y))
fprintf('     Proprietors'' income         %6.3f\n',mean(dopbuspinc./y))
fprintf('     Rental income               %6.3f\n',mean(dopbusrent./y))
fprintf('     Net interest                %6.3f\n',mean(dopbusnint./y))
fprintf('\n')
fprintf('   Government enterprises        %6.3f\n', ...
              mean((dgebus+dgebusdep)./y))
fprintf('     Capital consumption         %6.3f\n',mean(dgebusdep./y))
fprintf('     Compensation                %6.3f\n',mean(dgebuscomp./y))
fprintf('     Profits                     %6.3f\n',mean(dgebusgprof./y))
fprintf('\n')
fprintf('   Households                    %6.3f\n',mean((dhh+dhhdep)./y))
fprintf('     Capital consumption         %6.3f\n',mean(dhhdep./y))
fprintf('     IBT+transfers-subs.         %6.3f\n',mean(dhhibt./y))
fprintf('     Compensation                %6.3f\n',mean(dhhcomp./y))
fprintf('     Rental income               %6.3f\n',mean(dhhrent./y))
fprintf('     Net interest                %6.3f\n',mean(dhhnint./y))
fprintf('\n')
fprintf('   Nonprofits                    %6.3f\n',mean((dnp+dnpdep)./y))
fprintf('     Capital consumption         %6.3f\n',mean(dnpdep./y))
fprintf('     IBT+transfers-subs.         %6.3f\n',mean(dnpibt./y))
fprintf('     Compensation                %6.3f\n',mean(dnpcomp./y))
fprintf('     Rental income               %6.3f\n',mean(dnprent./y))
fprintf('     Net interest                %6.3f\n',mean(dnpnint./y))
fprintf('\n')
fprintf('   General Government            %6.3f\n',mean((dgg+dggdep)./y))
fprintf('     Capital consumption         %6.3f\n',mean(dggdep./y))
fprintf('     Compensation                %6.3f\n',mean(dggcomp./y))
fprintf('\n')
fprintf(' Statistical descrepancy         %6.3f\n',mean(desc./y))
fprintf('\n')
fprintf(' GDI                             %6.3f\n',mean(gdi./y))
fprintf('\n')
fprintf(' Net income, Rest of World       %6.3f\n',mean(row./y))
fprintf('   Compensation                  %6.3f\n',mean(rowcomp./y))
fprintf('   Profits                       %6.3f\n',mean(rowprof./y))
fprintf('   Net interest                  %6.3f\n',mean(rownint./y))
fprintf('\n')
fprintf('\n')
fprintf('\n')
fprintf(' DOMESTIC PRODUCT:\n')
fprintf('\n')
fprintf('   Private consumption           %6.3f\n',mean(pce./y))
fprintf('     Nondurable+Service          %6.3f\n',mean((pce-pced)./y))
fprintf('     Durable                     %6.3f\n',mean(pced./y))
fprintf('   Investment                    %6.3f\n',mean(gpdi./y))
fprintf('     Corporate tangible          %6.3f\n', ...
                mean((cgpdif+cgpdii)./y))
fprintf('       Fixed                     %6.3f\n',mean(cgpdif./y))
fprintf('       Inventory                 %6.3f\n',mean(cgpdii./y))
fprintf('     Noncorporate business       %6.3f\n', ...
                mean((ncgpdif+ncgpdii)./y))
fprintf('       Fixed                     %6.3f\n',mean(ncgpdif./y))
fprintf('       Inventory                 %6.3f\n',mean(ncgpdii./y))
fprintf('     Household                   %6.3f\n',mean(hhgpdif./y))
fprintf('     Nonprofits                  %6.3f\n',mean(npgpdif./y))
fprintf('   Government consumption        %6.3f\n',mean(govc./y))
fprintf('   Government investment         %6.3f\n',mean(govi./y))
fprintf('   Net Exports                   %6.3f\n',mean(netxm./y))
fprintf('\n')
fprintf(' GDP                             %6.3f\n',mean(gdp./y))
fprintf('\n')
fprintf('\n')
fprintf(' DOMESTIC CAPITAL:\n')
fprintf('\n')
fprintf('   Corporate                     %6.3f\n', ...
              mean((kcfix+kcinv+cland)./y))
fprintf('     Fixed                       %6.3f\n',mean(kcfix./y))
fprintf('     Inventory                   %6.3f\n',mean(kcinv./y))
fprintf('     Nonfarm land                %6.3f\n',mean(cland./y))
fprintf('   Noncorporate                  %6.3f\n', ...
              mean((kppfix+ktefix+kncinv+nfbland)./y))
fprintf('     Fixed,Proprietors           %6.3f\n',mean(kppfix./y))
fprintf('     Fixed,Tax-exempt            %6.3f\n',mean(ktefix./y))
fprintf('     Inventory                   %6.3f\n',mean(kncinv./y))
fprintf('     Nonfarm land                %6.3f\n',mean(nfbland./y))
fprintf('   Farm land (estimate)          %6.3f\n',.08)
fprintf('   Household                     %6.3f\n', ...
              mean((khhfix+hhland)./y))
fprintf('     Fixed                       %6.3f\n',mean(khhfix./y))
fprintf('     Land                        %6.3f\n',mean(hhland./y))
fprintf('   Nonprofits                    %6.3f\n',mean(knpfix./y))
fprintf('   Government                    %6.3f\n',mean(kgfix./y))
fprintf('   Consumer durables             %6.3f\n',mean(kcd./y))
fprintf('\n')
fprintf(' Total domestic capital          %6.3f\n',mean(k./y)+.08)
fprintf('\n')
fprintf('\n')
fprintf(' ADDENDUM:\n')
fprintf(' Consumer durable depreciation   %6.3f\n',mean(cddep./y))
fprintf(' IBT only                        %6.3f\n',mean(ibtonly./y))
fprintf('   Sales+other                   %6.3f\n',mean(stax./y))
fprintf('   Property                      %6.3f\n',mean(ptax./y))
fprintf(' Intermediate financial services %6.3f\n',mean(pbexp./y))
fprintf(' Net dividends                   %6.3f\n',mean(ndiv./y))      
fprintf(' Retained Earnings               %6.3f\n',mean(rete./y))
fprintf(' Dom. before-tax profits (no adj)%6.3f\n',mean(pbtd./y))
fprintf(' Dom. after-tax profits (no adj) %6.3f\n',mean(patd./y))
fprintf(' For. after-tax profits          %6.3f\n', ...
          mean((natprof-prof)./y))
fprintf(' Federal Reserve profits         %6.3f\n',mean(fedres./y))
fprintf(' Corporate tax liability         %6.3f\n',mean(ctax./y))
fprintf('\n')
fprintf('\n')
fprintf(' COMPARE NIPA and FA Fixed Investment\n')
fprintf('   Total \n')
fprintf('     NIPA                        %6.3f\n', ...
              mean((gpdif+pced+govi)./y))
fprintf('     FA                          %6.3f\n',mean((xfix+xcd)./y))
fprintf('   Corporate tangible\n')
fprintf('     NIPA                        %6.3f\n',mean(cgpdif./y))
fprintf('     FA                          %6.3f\n',mean(xcfix./y))
fprintf('   Noncorporate business\n')
fprintf('     NIPA                        %6.3f\n',mean(ncgpdif./y))
fprintf('     FA-Proprietors              %6.3f\n',mean(xppfix./y))
fprintf('     FA-Tax-exempt               %6.3f\n',mean(xtefix./y))
fprintf('   Government\n')
fprintf('     NIPA                        %6.3f\n',mean(govi./y))
fprintf('     FA                          %6.3f\n',mean(xgfix./y))
fprintf('   Consumer durables\n')
fprintf('     NIPA                        %6.3f\n',mean(pced./y))
fprintf('     FA                          %6.3f\n',mean(xcd./y))
fprintf('   Household\n')
fprintf('     NIPA                        %6.3f\n',mean(hhgpdif./y))
fprintf('     FA                          %6.3f\n',mean(xhhfix./y))
fprintf('   Nonprofits\n')
fprintf('     NIPA                        %6.3f\n',mean(npgpdif./y))
fprintf('     FA                          %6.3f\n',mean(xnpfix./y))
end;

%
% Save data files loaded into ./results.m
%
dtaup     = [t,taup*100];
usdiagnp  = [t,usdia./agnp*100];
fdiusgnp  = [t,fdius./agnp*100];
npignp    = [t,(rowprof+rownint-usdia+fdius)./agnp*100];
npiintgnp = [t,(rownint-dianint)./agnp*100];
cshare    = [t,(pce-pced+i*kcd+cddep-pbexp-cntax+govc+i*kgfix)./agdp*100];
nxgnp     = [t,netxm./agnp*100];
napfdi    = [t,usdia_cc./agnp*100];
napusdia  = [t,usdia_cc./agnp*100];
napfdius  = [t,fdius_cc./agnp*100];
napnet    = [t,(usdia_cc-fdius_cc)./agnp*100];
usdiaret  = [t,usdia./usdia_cc*100];
fdiusret  = [t,fdius./fdius_cc*100];

isave     = input('Update data input files? (1=yes, 0=no) ');
if isave;
save dtaup.dat dtaup -ascii
save dusdiagnp.dat usdiagnp -ascii
save dfdiusgnp.dat fdiusgnp -ascii
save dnpignp.dat npignp -ascii
save dnpiintgnp.dat npiintgnp -ascii
save dcshare.dat cshare -ascii
save dnxgnp.dat nxgnp -ascii
save dnapusdia.dat napusdia -ascii
save dnapfdius.dat napfdius -ascii
save dnapnet.dat napnet -ascii
save dusdiaret.dat usdiaret -ascii
save dfdiusret.dat fdiusret -ascii
end;
