Code covered by the BSD License

Toolkit on Econometrics and Economics Teaching

Hang Qian (view profile)

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))
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);

```