Here's an alternative PCA. It includes the 200 population averages that are the closest to the average coordinates of the new Amur River samples. Each population is connected with a line to its two closest neighbors.

Code:library(tidyverse) library(ggforce) # for geom_mark_hull t1=read.csv("https://drive.google.com/uc?export=download&id=1F2rKEVtu8nWSm7qFhxPU6UESQNsmA-sl",header=T,row.names=1,check.names=F) # ancient averages scaled t2=read.csv("https://drive.google.com/uc?export=download&id=1wZr-UOve0KUKo_Qbgeo27m-CQncZWb8y",header=T,row.names=1,check.names=F) # modern averages scaled t=rbind(t1,t2) amur=t1[grep("CHN_Amur_River.*BP",rownames(t1)),] mean=colMeans(amur) dist=apply(t,1,function(x)sqrt(sum((x-mean)^2))) t3=t[names(head(sort(dist),200)),] t=unique(rbind(amur,t3)) p=prcomp(t) pct=paste0(colnames(p$x)," (",sprintf("%.1f",p$sdev/sum(p$sdev)*100),"%)") p2=as.data.frame(p$x) k=cutree(hclust(dist(t)),16) p2$k=as.factor(k) set.seed(0) color=hcl(seq(15,375,length=nlevels(p2$k)+1)%>%head(-1),100,50)%>%sample xpc=1 ypc=2 xsym=sym(paste0("PC",xpc)) ysym=sym(paste0("PC",ypc)) dist=as.data.frame(as.matrix(dist(t))) seg=lapply(1:3,function(i)apply(dist,1,function(x)unlist(p2[names(sort(x)[i]),c(xpc,ypc)],use.names=F))%>%t%>%cbind(p2[,c(xpc,ypc)]))%>%do.call(rbind,.)%>%setNames(paste0("V",1:4)) ggplot(p2,aes(x=!!xsym,y=!!ysym))+ geom_segment(data=seg,aes(x=V1,y=V2,xend=V3,yend=V4),color="gray20",size=.1)+ geom_point(aes(color=k),size=.5)+ geom_mark_hull(aes(color=k,fill=k),concavity=1000,radius=unit(.15,"cm"),expand=unit(.15,"cm"),alpha=.2,size=.15)+ geom_text(aes(color=k),label=rownames(p2),size=2,vjust=-.7)+ # geom_text_repel(data=p2,aes(label=row.names(p2),color=as.factor(k)),size=2,max.overlaps=Inf,force=2,box.padding=.05,segment.size=.2,min.segment.length=.2)+ labs(x=pct[xpc],y=pct[ypc])+ coord_fixed()+ scale_x_continuous(breaks=seq(-1,1,.05))+ scale_y_continuous(breaks=seq(-1,1,.05))+ scale_color_manual(values=color)+ scale_fill_manual(values=color)+ theme( # aspect.ratio=1, axis.text.x=element_text(margin=margin(.2,0,0,0,"cm")), axis.text.y=element_text(angle=90,vjust=1,hjust=.5,margin=margin(0,.2,0,0,"cm")), axis.text=element_text(color="black",size=6), axis.ticks.length=unit(-.13,"cm"), axis.ticks=element_line(size=.3,color="gray80"), axis.title=element_text(color="black",size=8), legend.position="none", panel.background=element_rect(fill="white"), panel.border=element_rect(color="gray80",fill=NA,size=.6), panel.grid=element_blank() ) ggsave("a.png",width=12,height=12) system("mogrify -trim -bordercolor white -border 16 a.png")