Load the packages and data

library(vcdExtra)  # for mcaplot
library(ca)
data("PreSex", package="vcd")

Re-order the variables

We want the variables in the order G, P, E, M, with marital status last.

PreSex <- aperm(PreSex, 4:1)   # order variables G, P, E, M

multiple correspondence analysis, using the Burt matrix

presex.mca <- mjca(PreSex, lambda="Burt")
summary(presex.mca)
## 
## Principal inertias (eigenvalues):
## 
##  dim    value      %   cum%   scree plot               
##  1      0.149930  53.6  53.6  *************            
##  2      0.067201  24.0  77.6  ******                   
##  3      0.035396  12.6  90.2  ***                      
##  4      0.027365   9.8 100.0  **                       
##         -------- -----                                 
##  Total: 0.279892 100.0                                 
## 
## 
## Columns:
##                       name   mass  qlt  inr    k=1 cor ctr    k=2 cor ctr  
## 1 |             Gender:Men |   87  869  159 |  434 368 109 | -506 501 331 |
## 2 |           Gender:Women |  163  869   84 | -231 368  58 |  269 501 176 |
## 3 |       PremaritalSex:No |  192  763   61 | -251 714  81 |   66  49  12 |
## 4 |      PremaritalSex:Yes |   58  763  200 |  829 714 267 | -217  49  41 |
## 5 |     ExtramaritalSex:No |  221  681   29 | -149 596  33 |  -56  85  10 |
## 6 |    ExtramaritalSex:Yes |   29  681  221 | 1125 596 247 |  424  85  78 |
## 7 | MaritalStatus:Divorced |  119  794  128 |  368 450 108 |  322 344 184 |
## 8 |  MaritalStatus:Married |  131  794  117 | -336 450  98 | -294 344 168 |

the basic ca::plot() method doesn’t do a nice job of labeling

plot(presex.mca)

vcdExtra::mcaplot does a much nicer job

mcaplot(presex.mca)
# add a nice legend
cols <- c("blue", "red", "brown", "black")
legend("bottomright", legend=c("Gender", "PreSex", "ExtraSex", "Marital"), 
         title="Factor", 
       title.col="black",
       col=cols, text.col=cols, 
       pch=16:19, 
       bg="gray95", cex=1.2)

IycgLS0tDQojJyB0aXRsZTogIk11bHRpcGxlIENvcnJlc3BvbmRlbmNlIEFuYWx5c2lzOiBQcmUtIGFuZCBFeHRyYS1NYXJpdGFsIFNleCINCiMnIGF1dGhvcjogIk1pY2hhZWwgRnJpZW5kbHkiDQojJyBkYXRlOiAiYHIgZm9ybWF0KFN5cy5EYXRlKCkpYCINCiMnIG91dHB1dDoNCiMnICAgaHRtbF9kb2N1bWVudDoNCiMnICAgICB0aGVtZTogcmVhZGFibGUNCiMnICAgICBjb2RlX2Rvd25sb2FkOiB0cnVlDQojJyAtLS0NCg0KIysgZWNobz1GQUxTRQ0Ka25pdHI6Om9wdHNfY2h1bmskc2V0KA0KICB3YXJuaW5nID0gRkFMU0UsICAgIyBhdm9pZCB3YXJuaW5ncyBhbmQgbWVzc2FnZXMgaW4gdGhlIG91dHB1dA0KICBtZXNzYWdlID0gRkFMU0UNCikNCg0KIycgIyMgTG9hZCB0aGUgcGFja2FnZXMgYW5kIGRhdGENCg0KbGlicmFyeSh2Y2RFeHRyYSkgICMgZm9yIG1jYXBsb3QNCmxpYnJhcnkoY2EpDQpkYXRhKCJQcmVTZXgiLCBwYWNrYWdlPSJ2Y2QiKQ0KDQojJyAjIyBSZS1vcmRlciB0aGUgdmFyaWFibGVzDQojJyBXZSB3YW50IHRoZSB2YXJpYWJsZXMgaW4gdGhlIG9yZGVyIEcsIFAsIEUsIE0sIHdpdGggbWFyaXRhbCBzdGF0dXMgbGFzdC4NClByZVNleCA8LSBhcGVybShQcmVTZXgsIDQ6MSkgICAjIG9yZGVyIHZhcmlhYmxlcyBHLCBQLCBFLCBNDQoNCiMnICMjIG11bHRpcGxlIGNvcnJlc3BvbmRlbmNlIGFuYWx5c2lzLCB1c2luZyB0aGUgQnVydCBtYXRyaXgNCnByZXNleC5tY2EgPC0gbWpjYShQcmVTZXgsIGxhbWJkYT0iQnVydCIpDQpzdW1tYXJ5KHByZXNleC5tY2EpDQoNCiMnIHRoZSBiYXNpYyBjYTo6cGxvdCgpIG1ldGhvZCBkb2Vzbid0IGRvIGEgbmljZSBqb2Igb2YgbGFiZWxpbmcNCnBsb3QocHJlc2V4Lm1jYSkNCg0KIycgIyMgdmNkRXh0cmE6Om1jYXBsb3QgZG9lcyBhIG11Y2ggbmljZXIgam9iDQptY2FwbG90KHByZXNleC5tY2EpDQojIGFkZCBhIG5pY2UgbGVnZW5kDQpjb2xzIDwtIGMoImJsdWUiLCAicmVkIiwgImJyb3duIiwgImJsYWNrIikNCmxlZ2VuZCgiYm90dG9tcmlnaHQiLCBsZWdlbmQ9YygiR2VuZGVyIiwgIlByZVNleCIsICJFeHRyYVNleCIsICJNYXJpdGFsIiksIA0KCSAgICAgdGl0bGU9IkZhY3RvciIsIA0KICAgICAgIHRpdGxlLmNvbD0iYmxhY2siLA0KICAgICAgIGNvbD1jb2xzLCB0ZXh0LmNvbD1jb2xzLCANCiAgICAgICBwY2g9MTY6MTksIA0KICAgICAgIGJnPSJncmF5OTUiLCBjZXg9MS4yKQ0KDQo=