Introduction

This RMarkdown document performs a GMYC species delimitation analysis using the splits package in R.
GMYC (Generalized Mixed Yule-Coalescent) identifies species boundaries from a single-locus gene tree by detecting shifts in branching rates on an ultrametric and fully binary phylogenetic tree.

The method distinguishes between intraspecific (coalescent) and interspecific (Yule) diversification patterns.

Reference: T. Fujisawa’s GMYC website

Requirements

Before starting, ensure your input tree is:

Load Packages and Set Working Directory

# Set working directory - This is my working directory for this tutorial. Substitute with your own.
setwd("C:/1awinz/R_work/pseudochalceus")
getwd()
## [1] "C:/1awinz/R_work/pseudochalceus"
# Load necessary packages
library(ape)
library(splits)
## Loading required package: MASS
## Loading required package: paran

Load and Validate Tree

# Load the ultrametric tree - Change the name of the tree to what your tree is named.
tree <- read.nexus("1_MCC_tree_Aug_18_25.tree")

# Check if tree is suitable for GMYC
is.ultrametric(tree)  # Should be TRUE
## [1] TRUE
is.binary(tree)       # Should be TRUE
## [1] TRUE

Run GMYC Model

# Run GMYC model using the single threshold method
gmyc_result <- gmyc(tree, method = "single")
## node  T   loglik
## 2 -0.004327214 376.7427 
## 3 -0.003857657 376.9566 
## 4 -0.002630397 376.7507 
## 5 -0.002625162 376.6668 
## 6 -0.002009297 377.059 
## 7 -0.00141403 377.3936 
## 8 -0.001389607 377.3205 
## 9 -0.001276922 377.4747 
## 10 -0.001074432 377.8077 
## 11 -0.0009599549 378.6075 
## 12 -0.0007381195 378.9858 
## 13 -0.0007007941 378.9005 
## 14 -0.0005972165 378.7701 
## 15 -0.0005624536 379.4714 
## 16 -0.0004488043 378.6473 
## 17 -0.0004320271 377.9142 
## 18 -0.0004153024 377.3928 
## 19 -0.0003954619 376.3426 
## 20 -0.0003903773 376.1112 
## 21 -0.0003811927 375.4584 
## 22 -0.0003543288 376.4107 
## 23 -0.0003241173 376.5596 
## 24 -0.0003232253 376.04 
## 25 -0.0002960563 375.5184 
## 26 -0.0002846448 375.3241 
## 27 -0.0001549099 374.9398 
## 28 -0.0001509546 375.1129 
## 29 -0.0001497449 375.354 
## 30 -0.0001487261 375.6238 
## 31 -0.0001471876 375.888 
## 32 -0.0001466605 376.1222 
## 33 -0.0001461766 376.2989 
## 34 -0.0001424647 376.4079 
## 35 -0.0001423775 376.4788 
## 36 -0.0001409175 376.4887 
## 37 -0.0001404013 376.4628 
## 38 -0.0001396731 376.4076 
## 39 -0.0001395259 376.3399 
## 40 -0.0001390422 376.269 
## 41 -0.0001347165 376.2101 
## 42 -0.0001112986 376.2007 
## 
## Wed Aug 20 12:19:51 2025
## finish.
# Summary of the results
summary(gmyc_result)
## Result of GMYC species delimitation
## 
##  method: single
##  likelihood of null model:   376.3853
##  maximum likelihood of GMYC model:   379.4714
##  likelihood ratio:   6.172116
##  result of LR test:  0.04568167*
## 
##  number of ML clusters:  10
##  confidence interval:    8-12
## 
##  number of ML entities:  15
##  confidence interval:    9-17
## 
##  threshold time: -0.0005624536

Plot GMYC Output

# Plot GMYC output including threshold point and clusters
plot(gmyc_result)

Extract and Export Assignments

# Use spec.list to extract the sample-to-species mapping
assignments <- spec.list(gmyc_result)

# Display and optionally export the assignments
print(assignments)
##    GMYC_spec        sample_name
## 1          1             17_100
## 2          1             17_102
## 3          2             17_101
## 4          2              17_99
## 5          3             22_501
## 6          3             22_502
## 7          3             22_503
## 8          4               PL04
## 9          4               PL06
## 10         4               PL05
## 11         5    MECN_DP_4914_01
## 12         5    MECN_DP_4914_03
## 13         5    MECN_DP_4914_02
## 14         6    MECN_DP_4726_01
## 15         6    MECN_DP_4726_05
## 16         6    MECN_DP_4726_03
## 17         6    MECN_DP_4726_04
## 18         6    MECN_DP_4726_06
## 19         7    MECN_DP_4726_02
## 20         7    MECN_DP_4726_08
## 21         7    MECN_DP_4726_09
## 22         8 OR395303.1_Cube143
## 23         8  OR395318.1_Cube99
## 24         8 OR395306.1_Cube179
## 25         8 OR395313.1_Cube104
## 26         8 OR395308.1_Cube149
## 27         8  OR395319.1_Cube20
## 28         8  OR395314.1_Cube44
## 29         8   OR395316.1_Cube9
## 30         8 OR395311.1_Cube101
## 31         8  OR395315.1_Cube13
## 32         9 OR395305.1_Cube180
## 33         9  OR395309.1_Cube88
## 34         9  OR395317.1_Cube14
## 35         9 OR395307.1_Cube159
## 36         9 OR395312.1_Cube103
## 37        10 OR395304.1_Cube188
## 38        10 OR395310.1_Cube144
## 39        11             17_125
## 40        12            18V_723
## 41        13         Charsp_601
## 42        14    MECN_DP_4726_07
## 43        15               PL07
# Export to CSV
write.csv(assignments, "gmyc_spec_list.csv", row.names = FALSE)

Notes