Sequence analysis and descriptive statistics of clusters in R

Viewed 136

I am currently doing a sequence analysis, using the TraMineR package in R. However, i am having trouble finding out how to extract descriptive statistics for each cluster i get. Using the mvad dataset

mvad.seq <- seqdef(mvad, 17:86, alphabet = mvad.alphabet, states = mvad.scodes, 
    labels = mvad.labels, xtstep = 6)
clusterward1 <- agnes(dist.om1, diss = TRUE, method = "ward")
plot(clusterward1, which.plot = 2)
cl1.4 <- cutree(clusterward1, k = 4)
cl1.4fac <- factor(cl1.4, labels = paste("Type", 1:4))

How do I substract information about how many males are in each cluster, how many in each cluster are catholic, etc?

2 Answers

Although @Gilbert's answer already provides a very comprehensive solution, I had the impression that it might be somewhat intimidating for R novices who are only looking for a simple way of inspecting some cross tables.

If you really just want to obtain a simple cross table, you can use tablefunction suggested by Gilbert. It just requires two input vectors (e.g., mvad$male and cl1.4fac)

# Base R - cross table
table(cl1.4fac,
      mvad$male)

cl1.4fac  no yes
  Type 1 123 206
  Type 2 117  90
  Type 3  75  58
  Type 4  27  16

If you want to use the vector indicating cluster membership in additional analysis (e.g., regression), you might also want to consider adding it to the data frame containing the other relevant variables

# add cluster indicator to data
mvad$cluster <- cl1.4fac

If you are looking for a somewhat more appealing output, you can turn to additional libraries. In the example below I use {gtsummary}. Note that the following code is build on the premise that the cluster factor is part of mvad(see code above).

# nicely formatted cross table with gtsummary
library(gtsummary)

tbl_cross(
  data = mvad,
  row = cluster,
  col = male,
  percent = "row"
)

Cross table created with tbl_cross

With table you get the distribution of categorical variable. Applying table on each cluster you get the distribution by cluster. I illustrate below with a full working example for two variables (male and catholic).

library(TraMineR)
library(cluster)
data(mvad)
mvad.alphabet <- c("employment", "FE", "HE", "joblessness", "school", 
                   "training")
mvad.labels <- c("employment", "further education", "higher education", 
                 "joblessness", "school", "training")
mvad.scodes <- c("EM", "FE", "HE", "JL", "SC", "TR")
mvad.seq <- seqdef(mvad, 17:86, alphabet = mvad.alphabet, states = mvad.scodes, 
                   labels = mvad.labels, xtstep = 6)

dist.om1 <- seqdist(mvad.seq, method="OM", sm="INDELSLOG")
clusterward1 <- agnes(dist.om1, diss = TRUE, method = "ward")
cl1.4 <- cutree(clusterward1, k = 4)
cl1.4fac <- factor(cl1.4, labels = paste("Type", 1:4))

res <- list()
for (i in levels(cl1.4fac)){
  res[[i]] <- list()
  for (j in c("male","catholic")){
    res[[i]][[j]] <- table(mvad[cl1.4fac==i,j])
  }
}

Here is what you obtain

res

# $`Type 1`
# $`Type 1`$male
#
#  no yes
# 123 206
#
# $`Type 1`$catholic
#
#  no yes
# 183 146
#
#
# $`Type 2`
# $`Type 2`$male
#
#  no yes
# 117  90
#
# $`Type 2`$catholic
#
#  no yes
# 105 102
#
#
# $`Type 3`
# $`Type 3`$male
#
#  no yes
#  75  58
#
# $`Type 3`$catholic
#
#  no yes
#  66  67
#
#
# $`Type 4`
# $`Type 4`$male
#
#  no yes
#  27  16
#
# $`Type 4`$catholic
#
#  no yes
#  14  29 

res is a list of lists. You get for example the number of women in cluster 2 with

res[["Type 2"]][["male"]][1]
#  no 
# 117
Related