Figure legends:

FIGURE 1 Sampling sites and geographic distribution of the 12 ribotypes (A) and 19 chlorotypes (B) found in 10 populations ofPrimulina heterotricha from Hainan Island.

FIGURE 2 (A) Molecular phylogeny and (B) network of the 12 ribotypes from Primulina heterotricha based on nrITS data, numbers on the branches showed the Posterior probability (PP)/Bootstrap value (BS). The relative sizes of the circles in the network are proportional to the ribotype frequencies. ▲, vicariance events detected by S-DIVA. The dashed lines represent the isolation of two rivers in Hainan Island, and the size of dashed lines was correspondent to the magnitude of isolation. The red frames indicate different clades.

FIGURE 3 (A) Chronogram for 19 chlorotypes from Primulina heterotricha and outgroup based on combined trn L-trn F andycf 1b data, analysed using a relaxed molecular clock approach, posterior probabilities (PP) are showed above branches. (B) Network of 19 chlorotypes from Primulina heterotricha based on combinedtrn L-trn F and ycf 1b data, the relative sizes of the circles in the network are proportional to the chlorotypes frequencies, and missing chlorotypes are represented by small black spot. ▲, vicariance events detected by S-DIVA. The dashed lines represent the isolation of two rivers in Hainan Island, and the size of dashed lines was correspondent to the magnitude of isolation. The red frames indicate different clades.

FIGURE 4 Histogram of the STRUCTURE for all individuals at K=3 based on ITS1/2, the size of arrow lines were correspondent to the isolation level of rivers.

FIGURE 5 Results of Barrier analysis based on the (A) ITS1/2 and (B) combined trn L-trn F and ycf 1b, indicating the spatial separation of Primulina heterotricha populations, Delaunay triangulation (black lines) and inferred barriers (red, green and blue lines) separating the different original regions which showing the geographic location of the genetic barrier, black spots represent corresponding population location. NW, northwest of Changhua R., including populations BW, YJ, EX, WX and YG; SE, southeast of Changhua R., including populations WZ, QX, XA and QL; SW, southwest of Changhua R., including population JF.

FIGURE S1 The results graph of the relationship between genetic and geographic distance for 10 populations based on the (A) nrITS and (B) cpDNA.

FIGURE S2 Bayesian inference analysis of nrITS for the determination of the most likely number of genetic clusters (K).

FIGURE S3 Graphical output from S-DIVA. (A) Graphical results of ancestral distributions at each node of the phylogeny of chlorotypes from P. heterotricha obtained by S-DIVA. Alternative ancestral ranges of nodes I and II (with frequency of occurrence) are shown in pie chart form. Posterior probabilities (PP) are indicated by the sides of pie charts. (B) Colour blocks to possible ancestral ranges at different nodes. NW, northwest of Changhua R., including populations BW, YJ, EX, WX and YG; SE, southeast of Changhua R., including populations WZ, QX, XA and QL; SW, southwest of Changhua R., including population JF. ▲, possible vicariance events detected by S-DIVA.

FIGURE S4 Mismatch distribution analysis plots across all populations of (A) nrITS and (E) cpDNA and clade NW region of (B) nrITS and (F) cpDNA, clade SE region of (C) nrITS and (G) cpDNA, clade SW region of (D) nrITS and (H) cpDNA.