Brief history

Wheeler (1995) introduced a matrix plot to show how clade recovery varied across analytical parameters. Each cell represented a combination of gap cost and transversion-transition ratio, and was marked according to whether the clade was recovered as monophyletic, left unresolved, or resolved as nonmonophyletic under that parameter set. By plotting these matrices for individual clades across the full parameter space, Wheeler showed which groups were robust to analytical choices and which were sensitive to them. Giribet (2003) named these plots “Navajo rugs” and drew a distinction between nodal support and nodal stability. Sanders (2010) automated the approach with Cladescan (Perl), and Machado (2015) extended it with YBYRÁ (Python). Both produce standalone SVG files per node that have to be placed onto the tree by hand in a vector editor, and neither is maintained today. Both were also built for sensitivity analysis within a single analytical framework, varying parameters and cost schemes, rather than comparing trees from different inference pipelines that report support in incompatible formats (UFBoot2, SH-aLRT, posterior probability, ASTRAL LPP). phylorug moves this idea into R: it reads trees from any pipeline, bins support values against thresholds set for their own metric rather than forcing everything onto one scale, and draws the rug plot on the reference tree directly, with no manual placement and no outside software.

Load the data

sample_trees is a named list of five phylo objects, all sharing the same 15 taxa. Pick one as the backbone and use the rest as comparisons. In the original study 70p_uce is used as a backbone tree.

names(sample_trees)
#> [1] "70p_ASTRAL_partition_entropy" "70p_ASTRAL_uce"              
#> [3] "70p_ghost"                    "70p_partition_entropy"       
#> [5] "70p_uce"

backbone <- sample_trees[["70p_uce"]]
others   <- sample_trees[names(sample_trees) != "70p_uce"]

Check taxon consistency

check_taxa() compares the tip labels of the backbone against every comparison tree and reports any mismatches, missing taxa, extra taxa, or spelling differences. Running it before node_presence_matrix() is good practice, especially the first time you work with a new dataset. That said, node_presence_matrix() will stop with a clear error if a comparison tree is missing a backbone taxon, so if you already know your trees share the same taxa you can skip straight to building the matrix :

check_taxa(backbone, others)
#> All 4 comparison trees share the same 15 taxa as the backbone.
#> [1] TRUE
#> attr(,"diagnostics")
#>                     comparison    status n_taxa missing extra
#> 1 70p_ASTRAL_partition_entropy identical     15              
#> 2               70p_ASTRAL_uce identical     15              
#> 3                    70p_ghost identical     15              
#> 4        70p_partition_entropy identical     15

All trees share the same 15 taxa, so we can proceed. If a mismatch is reported the helper function prune_to_shared() can be used to get identical taxa set.

Build the node presence matrix

node_presence_matrix() is the core function of the pipeline. It walks every internal node of the backbone, checks whether each comparison tree recovered that clade, and records its support value.

It returns a named list: a presence matrix (1 = clade recovered, 0 = absent), plus one or more support matrices holding the support values pulled from each tree’s node labels. By default it reads the first value in each label. If your trees carry compound labels like 100/98 (SH-aLRT/UFBoot2 from IQ-TREE, for example), support_col picks which one to use: support_col = 1 for the first, support_col = 2 for the second, or support_col = c(1, 2) to pull both at once and makes two support matrices.

Two more arguments belong here. support_type is a named character vector telling node_presence_matrix() which support metric each comparison tree uses ("ufboot", "lpp", "sh_alrt", "posterior"). It gets stored with the result and plot_phylorug() reads it automatically in support mode, so you only declare it once, here.

pool_threshold matters when a comparison tree file holds several equally optimal trees (a pool) instead of one. It sets how many of them need to recover a clade before it counts as present: 1.0 (the default) a strict consensus; 0.5 needs just a majority; 0 skips the cutoff entirely and records the raw fraction that recovered it. sample_trees is all single phylo objects, no pools, so pool_threshold doesn’t come into play here.

We declare support_type now so plot_phylorug() can bin each tree’s support values against the right thresholds for its own metric:

support_type <- c(
  "70p_ASTRAL_partition_entropy" = "lpp",
  "70p_ASTRAL_uce"               = "lpp",
  "70p_ghost"                    = "ufboot",
  "70p_partition_entropy"        = "ufboot"
)

npm <- node_presence_matrix(backbone, others, support_col = c(1, 2), support_type = support_type)

The result is a named list. npm$presence is a matrix of 1s and 0s. npm$support_1 and npm$support_2 hold the raw support values. The list also carries support_type as an attribute, and all of these are what plot_phylorug() uses to draw the rug: presence determines whether a cell is filled or empty; support determines its shade.

Plot the rug

With the node presence matrix ready, plot_phylorug() draws the rug on the backbone. In its simplest form, all arguments are defaults: mode = "presence", dot_identical = TRUE (unanimous nodes get a black dot), legend = TRUE, show_support = TRUE (the backbone’s own support values are shown alongside the rug), and rug_position = "inside". The internal scaling engine sizes the canvas, fonts, and cell dimensions automatically based on tree size. Here show_support = FALSE keeps this first figure to presence and absence only.

plot_phylorug(backbone, npm, show_support = FALSE)

For publication-quality output, setting file, width, and height is strongly recommended. All other arguments are worth exploring to fine-tune your figure, but these three matter most.

plot_phylorug(backbone, npm, file = "presence.pdf", width = 12, height = 8)

Support mode

Now switch to support mode to see how strongly each tree supports each clade. Cells are shaded by binned support strength, darker means stronger. Support values are binned against thresholds specific to their own metric: a tree with LPP 0.95 and a tree with LPP 0.65 are compared to LPP thresholds, while a tree carrying UFBoot2 values is binned against UFBoot2 thresholds independently. phylorug never cross-compares values from different support metrics, because LPP 0.95 and UFBoot2 95 do not measure the same thing despite looking numerically similar.

Since we already declared support_type when building the node presence matrix, plot_phylorug() reads it automatically from the stored attribute and applies the correct metric-specific thresholds.

If support_type is not declared, plot_phylorug() falls back to a universal threshold scale: values greater than 1 are divided by 100 to bring them onto a 0–1 scale, and all trees are binned against the same set of thresholds regardless of metric. This is a reasonable starting point for quick exploration, but declaring support_type is always preferred because it ensures each value is evaluated against the thresholds published for its own method.

plot_phylorug(
  backbone, npm,
  mode         = "support")

In support mode a dot carries more information than in presence mode. It marks a clade that every comparison analysis recovers and rates very-high support against its own metric’s thresholds. A clade recovered by every analysis but weaker in even one, or with no usable support value, is drawn as a rug plot instead, so a dot never hides uneven support. Because the rule uses the same thresholds that shade the cells, changing thresholds also changes which nodes earn a dot (see Customisation).

Read the trees

read_trees() scans a directory for tree files and returns a named list:

tree_dir <- system.file("extdata", "culicomorpha", package = "phylorug")
trees    <- read_trees(tree_dir)
#> Read 10 analyses (10 trees) from: /tmp/RtmptPn09G/Rinst9ee7ceae29b/phylorug/extdata/culicomorpha
names(trees)
#>  [1] "Matrix1_LG_C20_F_R"           "Matrix1-kpi_ASTRAL"          
#>  [3] "Matrix1-kpi_partitioning"     "Matrix1-kpi_PMSF(H1.guide)"  
#>  [5] "Matrix1-kpi_PMSF(H2.guide)"   "Matrix1-smart_ASTRAL"        
#>  [7] "Matrix1-smart_partitioning"   "Matrix1-smart_PMSF(H1.guide)"
#>  [9] "Matrix1-smart_PMSF(H2.guide)" "Matrix2_partitioning"

Root and remove the outgroups

The Culicomorpha trees include three outgroup taxa from outside the infraorder: Phlebotomus chinensis and Clogmia albipunctata (Psychodidae) and Coboldia fuscipes (Scatopsidae). Root on all three, then drop them. The order matters, you must root while the outgroup tips are still present, because ape::root() needs to find them in the tree:

trees <- lapply(trees, function(tr) {
  tr <- ape::root(tr, outgroup = c("Coboldia_fuscipes",
                                   "Phlebotomus_chinensis",
                                   "Clogmia_albipunctata"),
                  resolve.root = TRUE)
  ape::drop.tip(tr, c("Coboldia_fuscipes",
                      "Phlebotomus_chinensis",
                      "Clogmia_albipunctata"))
})

The Culicomorpha trees already carry species names as tip labels, so no translation is needed here. For trees with specimen codes, see translate_tips() in the beetles example below.

Choose a backbone

Pick one analysis as the backbone. The remaining trees are comparisons:

backbone <- trees[["Matrix1-kpi_PMSF(H1.guide)"]]
others   <- trees[names(trees) != "Matrix1-kpi_PMSF(H1.guide)"]

Check taxon consistency

check_taxa(backbone, others)
#> All 9 comparison trees share the same 43 taxa as the backbone.
#> [1] TRUE
#> attr(,"diagnostics")
#>                     comparison    status n_taxa missing extra
#> 1           Matrix1_LG_C20_F_R identical     43              
#> 2           Matrix1-kpi_ASTRAL identical     43              
#> 3     Matrix1-kpi_partitioning identical     43              
#> 4   Matrix1-kpi_PMSF(H2.guide) identical     43              
#> 5         Matrix1-smart_ASTRAL identical     43              
#> 6   Matrix1-smart_partitioning identical     43              
#> 7 Matrix1-smart_PMSF(H1.guide) identical     43              
#> 8 Matrix1-smart_PMSF(H2.guide) identical     43              
#> 9         Matrix2_partitioning identical     43

Build the matrix and plot

support_type <- c(
  "Matrix1-kpi_ASTRAL"            = "lpp",
  "Matrix1-smart_ASTRAL"          = "lpp",
  "Matrix1-kpi_partitioning"      = "sh_alrt",
  "Matrix1-kpi_PMSF(H2.guide)"    = "sh_alrt",
  "Matrix1_LG_C20_F_R"            = "sh_alrt",
  "Matrix1-smart_partitioning"    = "sh_alrt",
  "Matrix1-smart_PMSF(H1.guide)"  = "sh_alrt",
  "Matrix1-smart_PMSF(H2.guide)"  = "sh_alrt",
  "Matrix2_partitioning"          = "sh_alrt"
)

npm <- node_presence_matrix(backbone, others, support_type = support_type)

Presence mode

plot_phylorug(backbone, npm, show_support = FALSE)

With 9 comparison trees, phylorug automatically arranges the rug grid (3 columns by 3 rows in this case). The position legend maps each cell position to an analysis name.

Support mode

Here we use show_support = TRUE to overlay the backbone’s own node labels in red for cross-referencing, with show_support_idx = 1 to display only the first support metric instead of the full compound label (e.g. showing 100 rather than 100/100). cell_scale = 0.3 shrinks the rug cells slightly for a denser tree. All-white rugs are hidden by default, so only contested and stable nodes remain visible.

plot_phylorug(backbone, npm,
              mode               = "support",
              show_support       = TRUE,
              show_support_idx   = 1,
              support_label_cex  = 0.35,
              support_label_col  = "red",
              cell_scale         = 0.3,
              rug_position       = "outside")

What the rug reveals

The rug immediately separates stable from contested regions of the Culicomorpha tree. Most nodes carry black dots, meaning every comparison analysis recovered the clade and rated it very-high support. All family-level groupings such as Culicidae, Chironomidae, Ceratopogonidae, and the Simuliidae + Thaumaleidae clade, show complete agreement across every inference method and dataset.

The main exception is the node uniting Chironomidae and Ceratopogonidae as sister groups, the most contentious relationship in culicomorph systematics. Eight of nine comparison analyses recover this clade with strong support, but the site-homogeneous partitioned analysis of the kpi-trimmed matrix (Matrix1-kpi_partitioning) fails to recover it entirely. This is consistent with the finding of Fu et al. (2025) that site-homogeneous models are more prone to missing this relationship, and the pattern is immediately visible in the rug without consulting multiple tree files.

Two additional areas show disagreement. Near the base of Chironomidae, the rug plots at the node uniting Podonomus and Parochlus steinenii show mixed shading, indicating that the internal arrangement of basal chironomid lineages varies across methods. Near the tips, the node at Nilodorum + Dicrotendipes (backbone support 75.0) shows grey shading, moderate support in the backbone that is not unanimously endorsed by all comparison analyses.

A few cells appear red at nodes within Chironomidae. Red means the analysis recovered the clade but no support value could be parsed from the node label. In this dataset, the original tree files do not carry support values at those nodes, so phylorug correctly marks them as not computed.

In contrast, the placement of Dixidae relative to the core Culicoidea clade is stable across all methods, a finding that would require checking each tree individually without phylorug.

This is the core value of phylorug: patterns that required opening 10 separate tree files side by side in the original study are condensed onto a single figure. A reader can immediately see which nodes are robust to analytical choice and which deserve further investigation.

Read the trees

tree_dir <- system.file("extdata", "beetles_70p", package = "phylorug")
trees    <- read_trees(tree_dir)
#> Read 5 analyses (5 trees) from: /tmp/RtmptPn09G/Rinst9ee7ceae29b/phylorug/extdata/beetles_70p
names(trees)
#> [1] "70p_ASTRAL_partition_entropy" "70p_ASTRAL_uce"              
#> [3] "70p_ghost"                    "70p_partition_entropy"       
#> [5] "70p_uce"

Root and drop outgroups

The beetle trees include two outgroup taxa, NicorbUCE and NicvesUCE. We will root on them first, then drop them. As before, rooting must happen while the outgroup tips are still in the tree.

trees <- lapply(trees, function(tr) {
  tr <- ape::root(tr, outgroup = c("NicorbUCE", "NicvesUCE"),
                  resolve.root = TRUE)
  ape::drop.tip(tr, c("NicorbUCE", "NicvesUCE"))
})

Translate tip labels

The trees still carry specimen codes at this point. A CSV file biogeo.csv included with the package maps each code to its species name. Let’s load it and check the first few rows:

biogeo_path <- system.file("extdata", "beetles_70p", "biogeo.csv",
                           package = "phylorug")
biogeo <- read.csv(biogeo_path)
head(biogeo)
#>   specimen_code                                    species_name
#> 1   STL10208208                    Amietina_larrochei__STL10208
#> 2     AmppriUCE              Amphistomus_primonactus_Amppri_UCE
#> 3    STL1009595                   Anonychonitis_freyi__STL10095
#> 4   Anop1COL892                         Anoplostethus_sp_COL892
#> 5    ApimmST003                         Aphodius_immundus_ST003
#> 6   STL10140140 Apotolamprus_aff_ambohitsitondronensi__STL10140

Now let’s translate. The ordering rule from earlier applies here. First root and drop outgroups before translating, because translate_tips() replaces the original codes and ape::root() would no longer find ” NicorbUCE ” after translation.

trees <- translate_tips(trees, biogeo,
                        from_col = "specimen_code",
                        to_col   = "species_name")
#> 70p_ASTRAL_partition_entropy: 43 tips translated, 0 unchanged
#> 70p_ASTRAL_uce: 43 tips translated, 0 unchanged
#> 70p_ghost: 43 tips translated, 0 unchanged
#> 70p_partition_entropy: 43 tips translated, 0 unchanged
#> 70p_uce: 43 tips translated, 0 unchanged

Choose a backbone and check taxa

backbone <- trees[["70p_uce"]]
others   <- trees[names(trees) != "70p_uce"]
check_taxa(backbone, others)
#> All 4 comparison trees share the same 43 taxa as the backbone.
#> [1] TRUE
#> attr(,"diagnostics")
#>                     comparison    status n_taxa missing extra
#> 1 70p_ASTRAL_partition_entropy identical     43              
#> 2               70p_ASTRAL_uce identical     43              
#> 3                    70p_ghost identical     43              
#> 4        70p_partition_entropy identical     43

Build the node presence matrix

The beetle trees carry compound node labels (SH-aLRT/UFBoot2 for IQ-TREE trees). Extract both support values with support_col = c(1, 2):

support_type <- c(
  "70p_partition_entropy"        = "ufboot",
  "70p_ghost"                    = "ufboot",
  "70p_ASTRAL_uce"               = "lpp",
  "70p_ASTRAL_partition_entropy" = "lpp"
)

npm <- node_presence_matrix(backbone, others, support_col = c(1, 2),
                            support_type = support_type)

Presence mode

At 43 taxa the tree is comfortably readable. All-white rugs are hidden by default, and writing to a file with explicit width and height is strongly recommended to avoid the distortion that GUI windows introduce.

plot_phylorug(backbone, npm)

Support mode

With 43 taxa, the canvas still fits a single page. The backbone tree carries compound IQ-TREE labels in the format SH-aLRT/UFBoot2 (e.g., 99.9/97). The show_support_idx argument controls which part of that label appears on the tree: show_support_idx = 1 shows only SH-aLRT, show_support_idx = 2 shows only UFBoot2, and a named vector like c("1" = "sh_alrt", "2" = "ufboot") shows both and names them in the legend. Here we use the named form so the legend identifies each metric:

# For publication output, add:
# file = "beetles_support.pdf"
# height = 12 
# width = 10
plot_phylorug(
  backbone,
  npm,
  mode               = "support",
  show_support       = TRUE,
  show_support_idx   = c("1" = "sh_alrt", "2" = "ufboot"),
  support_label_col  = "red",
  support_label_cex  = 0.38,
  cell_scale         = 0.35,
  cex                = 0.8
)

What the rug reveals

The 43-taxon subset corresponds exactly to the subtribe Onthophagina sensu lato as defined by Montanaro, Lopes et al. (2026, Fig. 2). Above we applied both presence mode and support mode to the same dataset. Presence mode shows which analyses recovered each clade; support mode goes further by showing how strongly each analysis backs each node.

Most nodes within Onthophagina carry black dots: all four comparison analyses recover the clade and rate it very-high support. The added value of support mode becomes visible at nodes where not all analyses recover the clade, or where they recover it with uneven support.

The clade containing species of Caccobius, Cleptocaccobius, and Onthophagus (backbone support 43.7/64 SH-aLRT/UFBoot2) is the most weakly supported node in the entire subtribe. No dot or rug appears at this node, meaning none of the comparison analyses recovered this clade. This confirms that the placement of Caccobius relative to Onthophagus is unresolved across methods.

Weak support and partial recovery are also visible deeper within the Onthophagus radiation. The node uniting O. giraffa + O. pilosus with O. bicavifrons + O. quadrimaculatus (backbone support 59.8/66) shows moderate-to-low support shading in the analyses that recover the clade, while others fail to recover it entirely. Within this group the sister pair O. quadrimaculatus + O. bicavifrons is recovered with strong support only by the backbone and the partitioned IQ-TREE analysis; the GHOST analysis supports it weakly, while both ASTRAL trees fail to recover the clade entirely. This is a clear coalescent-versus-concatenation disagreement that support mode makes visible. In presence mode alone, this node would look the same as any other partially recovered clade.

The Phalops group (Kurtops, Hamonthophagus, Digitonthophagus, Phalops, and O. probus) at the base of Onthophagina is recovered as a clade by all five analyses. However, the relationships inside the group tell a different story. The rug at the node uniting Hamonthophagus with Digitonthophagus + Phalops + O. probus (backbone support 63.8/73) shows mixed shading, meaning some analyses support this placement only weakly. One step deeper, the node grouping Digitonthophagus, Phalops, and O. probus (backbone support 99.9/97) also shows disagreement, high backbone support does not guarantee agreement across methods. Even the Digitonthophagus + Phalops pair (backbone support 100/100) is not recovered by the ASTRAL partitioned analysis. The Phalops group itself is stable, but the relationships among its members depend on the method.

These three cases: complete absence, partial weak recovery, and unanimous recovery with uncertain relationships within the group - capture different levels of conflict that support mode makes visible.

In support mode a black dot is a complete summary: the clade was recovered everywhere and every analysis has very-high support values. Nodes recovered everywhere but support is not uniformly very high are drawn as rug plots. To draw the grid at every node, including those that would otherwise get a dot, set dot_on_very_high = FALSE.

Dots in support mode

By default a dot marks only the nodes that every analysis recovers and rates very-high. To draw the full grid at every node instead, set dot_on_very_high = FALSE. This is the most conservative option, since no dot can stand in for detail. It affects support mode only; in presence mode, dot_identical controls dots.

plot_phylorug(backbone, npm, mode = "support", dot_on_very_high = FALSE)

Universal thresholds (no support_type)

If support_type is not declared when building the node presence matrix, plot_phylorug() falls back to universal thresholds. Values greater than 1 (such as UFBoot2 scores on the 0-100 scale) are divided by 100 internally so that all trees are binned against the same 0-1 cutoffs (0.95 / 0.80 / 0.50) regardless of metric. The legend displays these as percentages (>=95, 80-94, 50-79, <50) for readability:

npm_universal <- node_presence_matrix(backbone, others, support_col = c(1, 2))

plot_phylorug(
  backbone, npm_universal,
  mode = "support"
)
#> No `support_type` declared: values >1 are divided by 100 and binned against universal thresholds (0.95/0.80/0.50). To use metric-specific thresholds, pass `support_type` to `node_presence_matrix()`. To customise the universal thresholds, use the `thresholds` argument, e.g. `thresholds = list(universal = c(very_high = 0.95, high = 0.80, moderate = 0.50))`.

Grid dimensions

By default, phylorug chooses a roughly square grid. Override with n_rows and n_cols:

# rug with a flattened grid
plot_phylorug(backbone, npm, n_rows = 1, n_cols = 4)

Tree appearance

plot_phylorug() passes additional arguments through ... to ape::plot.phylo(). Common options:

plot_phylorug(backbone, npm, cex = 0.6, font = 3, edge.width = 1.5)

Including the backbone

By default, the backbone does not occupy a cell in the rug (it trivially recovers every clade, since it defines the topology). Set include_backbone = TRUE to add it as cell 1:

plot_phylorug(backbone, npm, include_backbone = TRUE)

In support mode, pass the backbone’s support metric directly to include_backbone so it gets binned against the correct thresholds:

plot_phylorug(backbone, npm,
              mode             = "support",
              include_backbone = "ufboot")

Custom thresholds

Override the default support bins with your own cutoffs:

my_thresholds <- list(
  ufboot = c(very_high = 98, high = 90, moderate = 70),
  lpp    = c(very_high = 0.99, high = 0.90, moderate = 0.70)
)

plot_phylorug(
  backbone, npm,
  mode       = "support",
  thresholds = my_thresholds
)

Each entry must name all three tiers (very_high, high, moderate); anything else stops early with an error naming the missing tier. The very_high cutoff also decides which nodes earn a dot in support mode: lowering UFBoot2’s from 95 to 90, for example, lets a node keep its dot when its only sub-threshold analysis scores 92.

Selecting which nodes to visualize

Sometimes you only want to look at a few clades, not the whole tree. nodes restricts the rug to just those, using the tip labels that define each clade:

plot_phylorug(
  backbone, npm,
  mode = "support",
  nodes = list(
    c("Gyronotus_pumilus_STL10003", "Scarabaeus_westwoodi_STL10034"),
    c("Nanos_dubitatus_STL5001", "Catharsius_sp._STL10033")
  )
)

Only these two clades get a dot or rug, every other node is left bare. If you already know the internal node ID, pass it directly as a numeric vector instead:

plot_phylorug(backbone, npm, nodes = 18)

But how do you find the right node ID in the first place? To find the nodes where analyses disagree, filter the presence matrix for rows with at least one 0, then look up the tips at each one:

disagree_ids <- as.integer(rownames(npm$presence))[
  apply(npm$presence, 1, function(row) any(row == 0))
]
disagree_ids
#> [1] 17 18 19 20 21 22 23

lapply(disagree_ids, function(id) ape::extract.clade(backbone, id)$tip.label)
#> [[1]]
#>  [1] "Gyronotus_pumilus_STL10003"        "Epactoides_hanski__STL29"         
#>  [3] "Nesosisyphus_pygmaeus__STL3"       "Sisyphus_muricatus__STL5"         
#>  [5] "Sisyphus_schaefferi_STL44"         "Nanos_dubitatus_STL5001"          
#>  [7] "Scarabaeus_westwoodi_STL10034"     "Kheper_nigroaeneus_STL10036"      
#>  [9] "Circellium_bacchus__STL10270"      "Epilissus_cuprarius_STL5011"      
#> [11] "Helictopleurus_fissicollis__STL28" "Catharsius_sp._STL10033"          
#> [13] "Copris_fidius_ST005"              
#> 
#> [[2]]
#> [1] "Gyronotus_pumilus_STL10003"        "Nanos_dubitatus_STL5001"          
#> [3] "Scarabaeus_westwoodi_STL10034"     "Kheper_nigroaeneus_STL10036"      
#> [5] "Circellium_bacchus__STL10270"      "Epilissus_cuprarius_STL5011"      
#> [7] "Helictopleurus_fissicollis__STL28" "Catharsius_sp._STL10033"          
#> [9] "Copris_fidius_ST005"              
#> 
#> [[3]]
#> [1] "Nanos_dubitatus_STL5001"           "Scarabaeus_westwoodi_STL10034"    
#> [3] "Kheper_nigroaeneus_STL10036"       "Circellium_bacchus__STL10270"     
#> [5] "Epilissus_cuprarius_STL5011"       "Helictopleurus_fissicollis__STL28"
#> [7] "Catharsius_sp._STL10033"           "Copris_fidius_ST005"              
#> 
#> [[4]]
#> [1] "Epilissus_cuprarius_STL5011"       "Helictopleurus_fissicollis__STL28"
#> [3] "Catharsius_sp._STL10033"           "Copris_fidius_ST005"              
#> 
#> [[5]]
#> [1] "Epilissus_cuprarius_STL5011"       "Helictopleurus_fissicollis__STL28"
#> [3] "Catharsius_sp._STL10033"          
#> 
#> [[6]]
#> [1] "Helictopleurus_fissicollis__STL28" "Catharsius_sp._STL10033"          
#> 
#> [[7]]
#> [1] "Nanos_dubitatus_STL5001"       "Scarabaeus_westwoodi_STL10034"
#> [3] "Kheper_nigroaeneus_STL10036"   "Circellium_bacchus__STL10270"

That shows you exactly which taxa sit at each disagreeing node, so you can pick the one you actually want to zoom in on with nodes, rather than guessing from a bare number.