% Data from
% http://www.pewglobal.org/files/2014/03/Morality-Topline-FINAL.pdf

%%%%%%%%%%%%%%%%%%%%%%%%%%%%
homo=[23 37 35 1 4 100
    30 15 50 1 4 100
    36 17 43 2 3 100
    36 14 50 0 0 100
    51 8 38 0 3 100
    42 19 31 2 5 100
    55 6 38 0 1 100
    24 45 24 6 2 100
    22 44 26 3 6 100
    56 14 27 1 3 100
    9 72 9 4 6 100
    4 78 12 3 3 100
    1 95 1 3 0 100
    2 95 3 0 0 100
    7 80 11 1 2 100
    1 94 4 0 1 100
    0 92 4 2 2 100
    27 43 25 1 4 100
    33 18 45 1 2 100
    13 61 17 3 6 100
    3 93 2 1 1 100
    38 31 25 2 5 100
    4 88 6 2 1 100
    1 85 3 0 11 100
    25 65 9 1 0 100
    18 57 21 3 1 100
    49 27 20 3 2 100
    32 51 9 2 5 100
    44 39 14 2 1 100
    45 32 13 6 5 100
    18 70 9 3 0 100
    36 40 15 4 5 100
    32 49 12 4 3 100
    1 98 1 1 0 100
    3 88 9 0 0 100
    1 85 11 1 1 100
    3 68 26 0 2 100
    18 62 12 5 3 100
    1 93 5 1 0 100];

contra=[52 7 36 1 4 100
    49 4 43 1 3 100
    62 3 32 2 1 100
    56 4 40 0 0 100
    70 1 28 0 1 100
    59 6 30 3 3 100
    64 2 32 1 1 100
    45 14 33 7 1 100
    47 17 24 9 3 100
    74 2 21 1 1 100
    62 7 19 5 7 100
    42 21 21 14 2 100
    35 7 36 21 1 100
    38 6 36 19 2 100
    34 13 31 20 3 100
    40 14 24 18 3 100
    27 39 17 15 2 100
    44 17 29 8 2 100
    55 4 39 1 1 100
    49 13 28 5 6 100
    65 10 15 10 1 100
    73 5 18 1 3 100
    21 40 23 14 2 100
    12 65 7 8 8 100
    60 29 9 1 1 100
    49 14 26 9 2 100
    69 6 21 2 2 100
    58 19 13 5 5 100
    79 8 11 1 1 100
    75 6 17 1 1 100
    63 23 11 2 1 100
    63 14 18 2 4 100
    82 5 8 3 1 100
    35 52 7 6 1 100
    47 33 11 7 2 100
    17 54 13 12 4 100
    51 28 11 6 3 100
    68 15 9 6 2 100
    51 38 6 5 1 100];

divorce=[33 22 36 6 2 100
    37 9 47 5 2 100
    41 9 42 6 2 100
    43 5 52 0 0 100
    58 7 32 2 1 100
    48 18 23 9 2 100
    57 4 37 2 0 100
    42 19 28 12 0 100
    35 22 24 16 2 100
    60 11 21 5 2 100
    46 22 13 16 3 100
    32 27 24 16 1 100
    51 7 30 11 1 100
    54 6 29 10 2 100
    48 9 28 13 2 100
    27 23 27 22 1 100
    16 46 22 15 1 100
    33 23 34 11 0 100
    43 9 43 3 1 100
    45 25 18 9 3 100
    20 57 10 13 0 100
    55 13 27 3 2 100
    16 45 22 15 1 100
    10 73 5 12 1 100
    25 67 7 1 0 100
    32 25 29 13 1 100
    60 17 20 2 1 100
    29 49 10 9 2 100
    62 24 10 4 0 100
    64 17 16 1 1 100
    30 48 9 11 1 100
    41 33 16 7 3 100
    58 24 8 9 1 100
    10 80 2 8 1 100
    17 59 13 10 1 100
    12 61 12 13 1 100
    35 32 16 17 1 100
    16 40 15 27 2 100
    9 76 6 10 0 100];


abort=[17 49 23 9 2 100
    26 26 37 7 4 100
    27 25 28 16 4 100
    38 14 47 0 0 100
    43 19 25 11 2 100
    25 41 11 20 3 100
    35 26 28 8 2 100
    16 54 11 18 1 100
    13 47 13 24 3 100
    49 18 18 12 3 100
    25 44 10 18 3 100
    18 52 14 15 1 100
    5 62 22 9 2 100
    4 53 36 5 3 100
    8 56 18 15 2 100
    7 54 13 24 1 100
    4 77 6 11 1 100
    25 35 23 15 1 100
    29 26 32 10 3 100
    29 37 20 9 4 100
    1 89 1 8 0 100
    44 28 19 6 2 100
    4 79 8 9 1 100
    1 85 4 4 6 100
    2 93 4 0 0 100
    19 58 12 10 1 100
    16 56 11 17 1 100
    6 83 4 7 1 100
    7 79 9 4 0 100
    15 64 9 8 3 100
    2 85 3 10 0 100
    16 63 12 8 2 100
    8 77 6 8 1 100
    2 92 2 3 1 100
    3 82 8 6 0 100
    2 80 10 7 1 100
    14 52 22 10 2 100
    10 61 8 19 2 100
    2 88 6 4 0 100];

alco=[32 16 46 4 1 100
    31 9 57 2 1 100
    38 9 47 5 2 100
    25 18 56 1 0 100
    41 14 40 4 2 100
    15 27 49 6 3 100
    27 22 46 4 1 100
    22 23 47 8 0 100
    37 17 34 10 2 100
    46 22 24 7 2 100
    26 44 13 15 2 100
    10 69 17 2 2 100
    5 75 13 5 3 100
    3 85 12 0 0 100
    25 51 21 1 2 100
    1 89 9 0 1 100
    5 85 6 3 1 100
    27 29 39 5 1 100
    36 10 46 7 1 100
    12 41 37 7 4 100
    2 85 8 4 0 100
    66 6 25 2 1 100
    10 70 15 4 1 100
    1 94 2 2 1 100
    34 39 24 2 0 100
    34 22 35 9 0 100
    26 42 27 3 1 100
    13 69 12 5 1 100
    29 47 20 3 0 100
    29 31 31 6 2 100
    8 82 8 2 0 100
    21 41 31 4 3 100
    36 39 17 7 1 100
    8 80 8 3 0 100
    9 59 24 7 0 100
    6 68 17 7 1 100
    9 63 24 3 1 100
    22 36 25 13 4 100
    19 57 14 8 1 100];

affair=[4 84 10 1 2 100
    4 76 17 1 2 100
    4 76 15 2 3 100
    12 47 40 0 0 100
    11 60 26 1 2 100
    10 64 16 7 3 100
    8 64 27 1 1 100
    5 79 12 5 0 100
    9 71 13 4 3 100
    17 66 10 4 3 100
    7 69 9 10 4 100
    2 94 2 1 1 100
    2 93 4 0 1 100
    1 93 5 0 1 100
    0 92 3 3 2 100
    1 94 4 1 1 100
    3 90 4 3 1 100
    11 73 13 3 1 100
    4 79 13 3 1 100
    8 74 12 3 3 100
    1 93 4 2 0 100
    12 69 14 1 4 100
    2 90 6 2 1 100
    0 92 3 0 4 100
    5 90 4 1 0 100
    8 81 8 2 0 100
    11 72 12 3 1 100
    7 84 6 3 1 100
    6 84 9 1 1 100
    13 67 16 3 2 100
    4 89 4 2 0 100
    9 73 12 4 2 100
    12 76 9 2 2 100
    10 86 2 2 0 100
    8 80 7 3 1 100
    7 77 12 3 1 100
    5 68 24 2 1 100
    11 65 11 11 2 100
    7 80 8 4 1 100];

gamble=[25 24 47 2 2 100
    19 23 55 2 2 100
    22 26 45 5 2 100
    31 13 56 0 0 100
    25 31 41 2 1 100
    12 33 48 4 2 100
    17 28 51 1 3 100
    7 52 37 4 1 100
    22 33 38 4 4 100
    11 65 19 2 2 100
    10 62 18 5 4 100
    5 80 12 1 1 100
    3 91 5 0 1 100
    2 91 6 0 1 100
    4 77 17 0 2 100
    1 91 7 0 1 100
    2 91 4 2 1 100
    16 51 30 2 1 100
    18 36 39 5 2 100
    8 55 28 5 4 100
    2 88 8 2 0 100
    38 31 28 2 2 100
    8 76 13 3 1 100
    0 95 3 0 1 100
    20 61 18 0 0 100
    7 78 12 2 1 100
    20 49 25 3 3 100
    11 61 17 4 7 100
    15 64 19 1 1 100
    15 48 27 6 3 100
    6 82 10 1 1 100
    13 49 30 4 4 100
    20 59 15 3 2 100
    3 91 3 3 0 100
    5 65 21 5 4 100
    1 77 16 2 3 100
    7 63 21 8 1 100
    18 40 27 10 5 100
    3 80 11 3 3 100];

unmarried=[29 30 36 1 3 100
    34 15 47 1 2 100
    38 13 44 3 2 100
    47 6 47 0 0 100
    57 6 34 1 2 100
    42 11 41 4 2 100
    52 8 39 0 0 100
    59 11 27 2 0 100
    44 22 21 9 4 100
    67 10 18 2 3 100
    40 30 10 15 4 100
    5 91 2 1 1 100
    2 90 6 1 0 100
    1 95 3 0 0 100
    5 81 8 3 2 100
    1 94 4 0 1 100
    4 89 3 2 1 100
    36 37 22 3 2 100
    34 15 47 1 2 100
    12 58 17 6 7 100
    1 97 1 1 0 100
    47 21 26 2 4 100
    3 87 5 3 1 100
    0 94 2 0 3 100
    20 71 7 1 0 100
    27 35 29 6 2 100
    45 22 27 3 2 100
    22 53 11 8 7 100
    47 35 14 3 0 100
    54 13 26 6 2 100
    22 56 14 7 0 100
    26 44 21 7 3 100
    61 21 13 4 2 100
    6 86 4 3 1 100
    4 79 13 3 1 100
    3 77 13 6 2 100
    6 63 26 4 1 100
    21 48 17 11 2 100
    9 77 9 5 1 100];

names={'United States',
    'Canada',
    'Britain',
    'France',
    'Germany',
    'Italy',
    'Spain',
    'Greece',
    'Poland',
    'Czech Republic',
    'Russia',
    'Turkey',
    'Egypt',
    'Jordan',
    'Lebanon',
    'Palest. ter.',
    'Tunisia',
    'Israel',
    'Australia',
    'China',
    'Indonesia',
    'Japan',
    'Malaysia',
    'Pakistan',
    'Philippines',
    'South Korea',
    'Argentina',
    'Bolivia',
    'Brazil',
    'Chile',
    'El Salvador',
    'Mexico',
    'Venezuela',
    'Ghana',
    'Kenya',
    'Nigeria',
    'Senegal',
    'South Africa',
    'Uganda'};

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%% Make triangular plots
clf

for issue=1:8
    subplot(3,3,issue)
    
switch(issue)
    case 1
        x=abort;
        tit='Abortion';
    case 2
        x=affair;
        tit='Extramarital affairs';
    case 3
        x=alco;
        tit='Alcohol'
    case 4
        x=contra;
        tit='Contraceptives';
    case 5
        x=divorce;
        tit='Divorce';
    case 6
        x=gamble;
        tit='Gambling';
    case 7
        x=homo;
        tit='Homosexuality';
    case 8
        x=unmarried
        tit='Unmarried sex';
end

z=zeros(size(x,1),1);
for j=1:size(x,1)
    a=x(j,1); b=x(j,2); c=x(j,3);
    z(j)=(1*a + exp(2*i*pi/3)*b + exp(2*i*2*pi/3)*c)/(a+b+c);
end
z=z*i;
plot(real(z), imag(z),'o');
hold on
%text(real(z), imag(z),char(names));
axis square
axis equal
axis off

% Border
plot(real(i*exp(i*2*(0:3)*pi/3)),imag(i*exp(i*2*(0:3)*pi/3)),'r-')
t={'Acceptable','Unacceptable','No issue'};
text(1.15*real(i*exp(i*2*(0:2)*pi/3)),1.15*imag(i*exp(i*2*(0:2)*pi/3)),char(t),'HorizontalAlignment','center','FontSize',6)

% Title
text(0,-.7,tit,'FontSize',10,'HorizontalAlignment','center')
end

pause

% Correlation matrix
clf
C=(corrcoef([abort(:,1),affair(:,1), alco(:,1), contra(:,1), divorce(:,1), gamble(:,1), homo(:,1), unmarried(:,1)]))
imagesc(C)
colorbar

 set(gca,'YTickLabel',{'Abortion','Affairs','Alcohol','Contra.','Divorce','Gambling','Homo.','Unmarried.'})
 set(gca,'XTickLabel',{'Abortion','Affairs','Alcohol','Contra.','Divorce','Gambling','Homo.','Unmarried.'})

pause

% PCA. 
% Uses just the first three variables of each issue, since the last two
% doesnt change much and makes it far harder to figure out the factors
clf
% Labels for variables
labels={'abort accept',
'abort unaccept',
'abort notissue',
'affair accept',
'affair unaccept',
'affair notissue',
'alco accept',
'alco unaccept',
'alco notissue',
'contra accept',
'contra unaccept',
'contra notissue',
'divorce accept',
'divorce unaccept',
'divorce notissue',
'gamble accept',
'gamble unaccept',
'gamle notissue',
'homo accept',
'homo unaccept',
'homo notissue',
'unmarried accept',
'unmarried unaccept',
'unmarried notissud',
};

A=[abort(:,1:3),affair(:,1:3), alco(:,1:3), contra(:,1:3), divorce(:,1:3), gamble(:,1:3), homo(:,1:3), unmarried(:,1:3)];
[COEFF, latent, explained] = pcacov(cov(A));

plot(A*COEFF(:,1),A*COEFF(:,2),'o')
text(A*COEFF(:,1),A*COEFF(:,2),char(names))

xlabel('Moral conservatism')
ylabel('Moral relevance')

pause

% 3D plot
clf
plot3(A*COEFF(:,1),A*COEFF(:,2),A*COEFF(:,3),'o')
text(A*COEFF(:,1),A*COEFF(:,2),A*COEFF(:,3),char(names))
xlabel('PCA1')
ylabel('PCA2')
zlabel('PCA3')
axis equal

pause
% Plot of factor loadings
pca=3; % factor to check
[dummy,ind]=sort(abs(COEFF(:,pca)))
plot(1:size(COEFF,1),COEFF(ind,pca))
text(1:size(COEFF,1),COEFF(ind,pca),char(labels(ind)))


pause
pca=3;

% Triangular plot with colors based on factor
clf
for issue=1:8
    subplot(3,3,issue)
    
switch(issue)
    case 1
        x=abort;
        tit='Abortion';
    case 2
        x=affair;
        tit='Extramarital affairs';
    case 3
        x=alco;
        tit='Alcohol'
    case 4
        x=contra;
        tit='Contraceptives';
    case 5
        x=divorce;
        tit='Divorce';
    case 6
        x=gamble;
        tit='Gambling';
    case 7
        x=homo;
        tit='Homosexuality';
    case 8
        x=unmarried
        tit='Unmarried sex';
end

z=zeros(size(x,1),1);
for j=1:size(x,1)
    a=x(j,1); b=x(j,2); c=x(j,3);
    z(j)=(1*a + exp(2*i*pi/3)*b + exp(2*i*2*pi/3)*c)/(a+b+c);
end
z=z*i;
scatter(real(z), imag(z),100,A*COEFF(:,pca),'filled');
hold on
%text(real(z), imag(z),char(names));
axis square
axis equal
axis off

plot(real(i*exp(i*2*(0:3)*pi/3)),imag(i*exp(i*2*(0:3)*pi/3)),'r-')
t={'Acceptable','Unacceptable','No issue'};
text(1.15*real(i*exp(i*2*(0:2)*pi/3)),1.15*imag(i*exp(i*2*(0:2)*pi/3)),char(t),'HorizontalAlignment','center','FontSize',6)

text(0,-.7,tit,'FontSize',10,'HorizontalAlignment','center')
end


pause
% Plot of PC2 and PC3 for high and low scorers of PC1
clf
hold on
flow=find(A*COEFF(:,1)<0); 
text(A(flow,:)*COEFF(:,2),A(flow,:)*COEFF(:,3),char(names(flow)),'FontSize',8,'Color',[1 1 1]*.7)
plot(A(flow,:)*COEFF(:,2),A(flow,:)*COEFF(:,3),'.','MarkerSize',30)

fhi=find(A*COEFF(:,1)>100); 
text(A(fhi,:)*COEFF(:,2),A(fhi,:)*COEFF(:,3),char(names(fhi)),'FontSize',8,'Color',[1 1 1]*.7)
plot(A(fhi,:)*COEFF(:,2),A(fhi,:)*COEFF(:,3),'r.','MarkerSize',30)

xlabel('PCA2')
ylabel('PCA3')