## R-scripts associated to the paper entitled : "Repeated morphological diversification in endemic Antarctic fishes of the genus Trematomus" published in the Belgian Journal of Zoology 152: 55-73. Frédérich et al. (2022).

######## Disparity analysis #########
## Ask for the used packages
library(geomorph);

### Read TPS file ########
TPS_data <- readland.tps("Head_Trematomus_ALL.TPS", specID="imageID");

### Read classifier ########
classifier <- read.table("classifier.txt", header=T, row.names=1);
Species<- classifier[,1];
Group1<- classifier[,2];
Group2<- classifier[,3];
Group3<- classifier[,4];
Group4<- classifier[,5];
Group5<- classifier[,6];
Group6<- classifier[,7];

###### GPA superimposition ########
GPanalyses <- gpagen(TPS_data,Proj=TRUE,PrinAxes=TRUE);

###### Calculate and Compare Level of disparity ########

data.frame<-geomorph.data.frame(GPanalyses, Group1 = Group1, Group2 = Group2, Group3 = Group3, Group4 = Group4, Group5 = Group5, Group6 = Group6); 

Disparity_level <- morphol.disparity(coords ~ log(Csize), groups = ~ GroupX, data=data.frame, iter = 9999); 
## "groups = ~ GroupX", replace GroupX by the grouping factor you want.  

Disparity_level

############################################################
############################################################
############################################################

######## Comparisons of covariance matrices #########
## Ask for the used packages
library(geomorph);
library(evolqg);

#### Generate VCV matrix ###
# Before the analyses, you need to build different TPS files where you grouped species.
# Hereafter, we provide the exemple of one comparison : "bernacchii/hansoni" Vs "eulepidotus/loennbergii/lepidorhinus".  

TPS_data <- readland.tps("Head_bernacchiiANDhansoni.tps",specID="imageID"); ##
GPanalyses <- gpagen(TPS_data,Proj=TRUE,PrinAxes=TRUE);
CV_matrix_group1<-GPanalyses$points.VCV

TPS_data2 <- readland.tps("Head_eulepidotus_loennbergii_lepidorhinus.tps", specID="imageID");
GPanalyses2 <- gpagen(TPS_data2,Proj=TRUE,PrinAxes=TRUE);
CV_matrix_group2<-GPanalyses2$points.VCV


test <- MatrixCompare(CV_matrix_group1, CV_matrix_group2);
test ### MatrixCompare provides results from various correlation methods. See descriptions of function in evolqg.