Code covered by the BSD License  

Highlights from
Toolkit on Econometrics and Economics Teaching

image thumbnail

Toolkit on Econometrics and Economics Teaching

by

 

Many MATLAB routines related to econometrics, statistics and introductory economics teaching.

PPF.m
% Written by Hang Qian, Iowa State University
% Contact me:  matlabist@gmail.com

clc
fid_write = fopen('Output.txt', 'wt');

%-------------------------
% ----Model setup------
%-------------------------
%  Two countries a, b
%  Two goods 1,2
%  Production costs: Pa1, Pa2, Pb1, Pb2
%  Domestic resource constraints Ra, Rb
Ra = str2double(get(H_Ra,'string'));
Pa1 = str2double(get(H_Pa1,'string'));
Pa2 = str2double(get(H_Pa2,'string'));

Rb = str2double(get(H_Rb,'string'));
Pb1 = str2double(get(H_Pb1,'string'));
Pb2 = str2double(get(H_Pb2,'string'));

Xa1_old = str2double(get(H_Xa1_old,'string'));
Xb1_old = str2double(get(H_Xb1_old,'string'));

Name_a = get(H_Name_a,'string');
Name_b = get(H_Name_b,'string');
Name_1 = get(H_Name_1,'string');
Name_2 = get(H_Name_2,'string');
close(H_figure)


%-------------------------
% -----AUTARKY ------
%-------------------------
% Domestic maximum output
Xa1_max = Ra / Pa1;
Xa2_max = Ra / Pa2;
Xb1_max = Rb / Pb1;
Xb2_max = Rb / Pb2;

% PPF in autarky
% Country a PPFXa1_aut  v.s.  Xa2_aut
% Country b PPFXb1_aut  v.s.  Xb2_aut
Xa1_aut = linspace(0,Xa1_max, Xa1_max+1)';
Xa2_aut = (Ra - Xa1_aut * Pa1) / Pa2;
Xb1_aut = linspace(0,Xb1_max, Xb1_max+1)';
Xb2_aut = (Rb - Xb1_aut * Pb1) / Pb2;

if Pa1 < Pb1
    disp([Name_a, ' has an absolute advantage in ', Name_1]);
    fprintf(fid_write,'%s has an absolute advantage in %s \n',Name_a, Name_1);
elseif Pa1 > Pb1
    disp([Name_b, ' has an absolute advantage in ', Name_1]);
    fprintf(fid_write,'%s has an absolute advantage in %s \n',Name_b, Name_1);
else
    disp(['Two countries has the same productivity in ', Name_1, ', hence no absolute advantage.'])
    fprintf(fid_write,'Two countries has the same productivity in %s hence no absolute advantage. \n', Name_1);
end

if Pa2 < Pb2
    disp([Name_a, ' has an absolute advantage in ', Name_2]);
    fprintf(fid_write,'%s has an absolute advantage in %s \n',Name_a, Name_2);
elseif Pa2 > Pb2
    disp([Name_b, ' has an absolute advantage in ', Name_2]);
    fprintf(fid_write,'%s has an absolute advantage in %s \n',Name_b, Name_2);
else
    disp(['Two countries has the same productivity in ', Name_2, ', hence no absolute advantage.'])
    fprintf(fid_write,'Two countries has the same productivity in %s hence no absolute advantage. \n', Name_2);
end
disp(' ')
fprintf(fid_write,'\n');

%-------------------------
% -----OPEN ECONOMY ------
%-------------------------

XT1_max = Xa1_max + Xb1_max;
ngrid = floor(XT1_max) + 1;
XT1 = linspace(0,XT1_max, ngrid)';

% Opportunity cost of producing Good 1: OCa, OCb
OCa = Pa1 / Pa2;
OCb = Pb1 / Pb2;
disp(['In ', Name_a,', the opportunity cost of producing one unit of ', Name_1, ' is  ', num2str(OCa), '  units of ', Name_2])
disp(['In ', Name_b,', the opportunity cost of producing one unit of ', Name_1, ' is  ', num2str(OCb), '  units of ', Name_2])
fprintf(fid_write,'In %s , the opportunity cost of producing one unit of  %s is %6.2f units of %s \n', Name_a, Name_1, OCa, Name_2);
fprintf(fid_write,'In %s , the opportunity cost of producing one unit of  %s is %6.2f units of %s \n', Name_b, Name_1, OCb, Name_2);
if OCa < OCb
    disp([Name_a, ' has a comparative advantage in ', Name_1, ' (', Name_b, ' has a comparative advantage in ', Name_2, ')'])
    fprintf(fid_write,'%s has a comparative advantage in %s  (%s has a comparative advantage in %s ) \n',Name_a, Name_1, Name_b, Name_2);
else
    disp([Name_b, ' has a comparative advantage in ', Name_1, ' (', Name_a, ' has a comparative advantage in ', Name_2, ')'])
    fprintf(fid_write,'%s has a comparative advantage in %s  (%s has a comparative advantage in %s ) \n',Name_b, Name_1, Name_a, Name_2);
end
disp(' ')
disp('In perfect specialization, therefore, ')
fprintf(fid_write,'\n In perfect specialization, therefore \n');

% Optimal allocation PPF:  XT1 v.s. XT2
if OCa < OCb
    Xa1 = min(XT1, repmat(Xa1_max, ngrid,1));
    Xa2 = (Ra - Xa1 * Pa1) / Pa2;
    Xb1 = XT1 - Xa1;
    Xb2 = (Rb - Xb1 * Pb1) / Pb2;
else
    Xb1 = min(XT1, repmat(Xb1_max, ngrid,1));
    Xb2 = (Rb - Xb1 * Pb1) / Pb2;
    Xa1 = XT1 - Xb1;
    Xa2 = (Ra - Xa1 * Pa1) / Pa2;
end
XT2 = Xa2 + Xb2;

% Perfect specialization
if OCa < OCb
    disp([Name_a, ' should specialize in ', Name_1, ', with production ', num2str(Xa1_max)])
    disp([Name_b, ' should specialize in ', Name_2, ', with production ', num2str(Xb2_max)])
    fprintf(fid_write,'%s should specialize in %s, with production %6.2f \n',Name_a, Name_1,Xa1_max);
    fprintf(fid_write,'%s should specialize in %s, with production %6.2f \n',Name_b, Name_2,Xb2_max);
else
    disp([Name_a, ' should specialize in ', Name_2, ', with production ', num2str(Xa2_max)])
    disp([Name_b, ' should specialize in ', Name_1, ', with production ', num2str(Xb1_max)])
    fprintf(fid_write,'%s should specialize in %s, with production %6.2f \n',Name_a, Name_2,Xa2_max);
    fprintf(fid_write,'%s should specialize in %s, with production %6.2f \n',Name_b, Name_1,Xb1_max);
end

% Plot PPF graph
figure(1)
plot(Xa1_aut, Xa2_aut, 'b', 'LineWidth',2);
title( ['PPF of ', Name_a ,' in Autarky'],  'fontsize',12);
xlim( [0, max(XT1)*1.1] );  ylim( [0, max(XT2)*1.1] );
xlabel(Name_1,'fontsize',12); ylabel(Name_2, 'fontsize',12);

figure(2)
plot(Xb1_aut, Xb2_aut,'r', 'LineWidth',2);
title( ['PPF of ', Name_b ,' in Autarky'], 'fontsize',12);
xlim( [0, max(XT1)*1.1] );  ylim( [0, max(XT2)*1.1] );
xlabel(Name_1, 'fontsize',12); ylabel(Name_2, 'fontsize',12);

figure(3)
plot(Xa1_aut, Xa2_aut,'b', 'LineWidth',2);
hold on
plot(Xb1_aut, Xb2_aut, 'r', 'LineWidth',2); 
plot(XT1, XT2, 'g', 'LineWidth',2);
title('PPF in the Open Economy', 'fontsize',12);
xlim( [0, max(XT1)*1.1] );  ylim( [0, max(XT2)*1.1] );
xlabel(Name_1, 'fontsize',12); ylabel(Name_2, 'fontsize',12);
legend(['PPF of ', Name_a], ['PPF of ', Name_b], 'PPF open economy', 'Location','NorthEast')
% Perfect specialization PPF
if OCa < OCb
    scatter(Xa1_max, Xb2_max)
    plot( linspace(0, Xa1_max, ngrid), Xb2_max, '--', 'LineWidth',0.5)
    plot( Xa1_max, linspace(0, Xb2_max, ngrid), '--', 'LineWidth',0.5)
    text(Xa1_max * 1.05, Xb2_max * 1.05, 'Perfect specialization')
else
    scatter(Xb1_max, Xa2_max)
    plot( linspace(0, Xb1_max, ngrid), Xa2_max, '--', 'LineWidth',0.5)
    plot( Xb1_max, linspace(0, Xa2_max, ngrid), '--', 'LineWidth',0.5)
    text(Xb1_max * 1.05, Xa2_max * 1.05, 'Perfect specialization')
end
hold off


%--------------------------------
disp(' ')
fprintf(fid_write,'\n');

% Autarky production bundle
if isempty(Xa1_old) || Xa1_old > Xa1_max
    disp('Warning: Autarky production not specified / not feasibe, use default value instead.')
    Xa1_old = 0.5 * Xa1_max;
end
Xa2_old = (Ra - Xa1_old * Pa1) / Pa2;
disp(['In autarky, ', Name_a, ' consumes ', num2str(Xa1_old), ' units of ', Name_1])
disp(['In autarky, ', Name_a, ' consumes ', num2str(Xa2_old), ' units of ', Name_2])
fprintf(fid_write,'In autarky, %s consumes %6.2f units of %s \n',Name_a, Xa1_old, Name_1);
fprintf(fid_write,'In autarky, %s consumes %6.2f units of %s \n\n',Name_a, Xa2_old, Name_2);
disp(' ')


if isempty(Xb1_old) || Xb1_old > Xb1_max
    disp('Warning: Autarky production not specified / not feasibe, use default value instead.')
    Xb1_old = 0.5 * Xb1_max;
end
Xb2_old = (Rb - Xb1_old * Pb1) / Pb2;
disp(['In autarky, ', Name_b, ' consumes ', num2str(Xb1_old), ' units of ', Name_1])
disp(['In autarky, ', Name_b, ' consumes ', num2str(Xb2_old), ' units of ', Name_2])
fprintf(fid_write,'In autarky, %s consumes %6.2f units of %s \n',Name_b, Xb1_old, Name_1);
fprintf(fid_write,'In autarky, %s consumes %6.2f units of %s \n\n',Name_b, Xb2_old, Name_2);
disp(' ')

XT1_old = Xa1_old + Xb1_old;
XT2_old = Xa2_old + Xb2_old;
index = (XT1 - XT1_old >= 0) & (XT2 - XT2_old >= 0);

% Mutual benefit region
figure(4)
plot(Xa1_aut, Xa2_aut,'b', 'LineWidth',2);
hold on
plot(Xb1_aut, Xb2_aut, 'r', 'LineWidth',2); 
plot(XT1, XT2, 'g', 'LineWidth',2);
scatter(XT1(index), XT2(index))
title('PPF After Trade', 'fontsize',12);
xlim( [0, max(XT1)*1.1] );  ylim( [0, max(XT2)*1.1] );
xlabel(Name_1, 'fontsize',12); ylabel(Name_2, 'fontsize',12);
legend(['PPF of ', Name_a], ['PPF of ', Name_b], 'PPF After Trade','Mutual benefit region',  'Location','NorthEast')
scatter(Xa1_old,Xa2_old);
scatter(Xb1_old,Xb2_old);

% Example: one mutual benefit plan
share = 0.5;
if OCa < OCb
    ind = (XT1 == Xa1_max) & index;
else
    ind = (XT1 == Xb1_max) & index;
end
if ~any(ind)
    ind = ceil(median(XT1(index)));
end
XT1_new = XT1(ind);
XT2_new = XT2(ind);
Consumption_a1 = Xa1_old + ceil( (XT1_new - XT1_old) * share);
Consumption_b1 = XT1_new - Consumption_a1;
Consumption_a2 = Xa2_old + ceil( (XT2_new - XT2_old) * share);
Consumption_b2 = XT2_new - Consumption_a2;
Production_a1 = Xa1(ind);
Production_b1 = Xb1(ind);
Production_a2 = Xa2(ind);
Production_b2 = Xb2(ind);
disp('One feasible mutual benefit plan: ')
disp(['For ', Name_1])
fprintf(fid_write,'One feasible mutual benefit plan: \nFor %s \n',Name_1);
if Production_a1 - Consumption_a1 > 0
   flag1  = 'exports';   flag2 = 'imports';
else
   flag1  = 'imports';   flag2 = 'exports'; 
end
disp( [Name_a, ' produces ', num2str(Production_a1), ' and ', flag1,' ', ...
    num2str(abs(Production_a1 - Consumption_a1)), ' , so it consumes ', num2str(Consumption_a1), ...
    ' , better than autarky consumption ', num2str(Xa1_old)] );
disp( [Name_b, ' produces ', num2str(Production_b1), ' and ', flag2,' ', ...
    num2str(abs(Production_b1 - Consumption_b1)), ' , so it consumes ', num2str(Consumption_b1), ...
    ' , better than autarky consumption ', num2str(Xb1_old)] );
fprintf(fid_write,'%s produces %6.2f and %s %6.2f , so it consumes %6.2f, better than autarky consumption %6.2f \n',...
    Name_a, Production_a1, flag1,abs(Production_a1 - Consumption_a1), Consumption_a1, Xa1_old);
fprintf(fid_write,'%s produces %6.2f and %s %6.2f , so it consumes %6.2f, better than autarky consumption %6.2f \n',...
    Name_b, Production_b1, flag2,abs(Production_b1 - Consumption_b1), Consumption_b1, Xb1_old);

disp(['For ', Name_2])
fprintf(fid_write,'For %s \n',Name_2);
disp( [Name_a, ' produces ', num2str(Production_a2), ' and ',flag2,' ',  ...
    num2str(abs(Production_a2 - Consumption_a2)), ' , so it consumes ', num2str(Consumption_a2), ...
    ' , better than autarky consumption ', num2str(Xa2_old)] );
disp( [Name_b, ' produces ', num2str(Production_b2), ' and ', flag1,' ', ...
    num2str(abs(Production_b2 - Consumption_b2)), ' , so it consumes ', num2str(Consumption_b2), ...
    ' , better than autarky consumption ', num2str(Xb2_old)] );
fprintf(fid_write,'%s produces %6.2f and %s %6.2f , so it consumes %6.2f, better than autarky consumption %6.2f \n',...
    Name_a, Production_a2, flag2,abs(Production_a2 - Consumption_a2), Consumption_a2, Xa2_old);
fprintf(fid_write,'%s produces %6.2f and %s %6.2f , so it consumes %6.2f, better than autarky consumption %6.2f \n',...
    Name_b, Production_b2, flag1,abs(Production_b2 - Consumption_b2), Consumption_b2, Xb2_old);
scatter(XT1_new, XT2_new)
text(XT1_new, XT2_new,'   One feasible plan')
hold off
disp(' ')
disp('Program finished. Thanks for watching trade PPF demonstration.')
fprintf(fid_write,'\nProgram finished. Thanks for watching trade PPF demonstration.');
fclose(fid_write);

Contact us