Auximetric analysis and the slow-fast life-history continuum in marine fish and the coelacanth
Bruno Ernande
11/05/2021
Setting the working directory.
Length-based auximetric analysis
Growth data
We start importing the data. These data were extracted from Fishbase using the package rfishbase on June 2020. Only species from brackish and salted waters were included. Von Bertalanffy growth parameters, namely the rate coefficient \(K\) and asymptotic length \(L_\infty\) were imported from the Popgrowth table. As recommended by Fishbase, stocks for which the estimated \(L_\infty\) was outside the range 0.66-1.33 times the maximum length ever observed \(L_{max}\) were considered dubious and thus excluded. Likewise, stocks whose estimates were outside auximetric ellipses were removed. Notice though that a large number of stocks were not tested and were kept in the analysis. Total length was used for all species except for Scombridaes for which fork length was considered.
growthdata<-read.csv("VB_Growth_data_Fishbase_2020_June.csv")
summary(growthdata)## AutoCtr SpecCode Species StockCode
## Min. : 139 Min. : 4 Length:5517 Min. : 4
## 1st Qu.: 2177 1st Qu.: 204 Class :character 1st Qu.: 408
## Median : 4703 Median : 1288 Mode :character Median : 1672
## Mean : 7815 Mean : 3828 Mean : 7248
## 3rd Qu.:10235 3rd Qu.: 4449 3rd Qu.: 4986
## Max. :23914 Max. :67314 Max. :58558
##
## Locality PopGrowthRef Sex Loo
## Length:5517 Min. : 9 Length:5517 Min. : 2.10
## Class :character 1st Qu.: 2081 Class :character 1st Qu.: 25.90
## Mode :character Median : 13304 Mode :character Median : 41.00
## Mean : 38231 Mean : 62.46
## 3rd Qu.: 74354 3rd Qu.: 77.20
## Max. :119692 Max. :1966.00
##
## TLinfinity K to Winfinity
## Min. : 2.10 Min. : 0.0070 Min. :-12.210 Min. : 0
## 1st Qu.: 26.80 1st Qu.: 0.1590 1st Qu.: -1.700 1st Qu.: 193
## Median : 41.97 Median : 0.2700 Median : -0.858 Median : 721
## Mean : 63.92 Mean : 0.4517 Mean : -1.200 Mean : 31059
## 3rd Qu.: 79.00 3rd Qu.: 0.5300 3rd Qu.: -0.300 3rd Qu.: 4526
## Max. :1966.00 Max. :12.9570 Max. : 1.100 Max. :20000000
## NA's :2883 NA's :3091
## Type Auxim LinfLmax LogKLogLoo
## Length:5517 Length:5517 Min. :0 Min. :-0.8340
## Class :character Class :character 1st Qu.:0 1st Qu.:-0.4215
## Mode :character Mode :character Median :0 Median :-0.3149
## Mean :0 Mean :-0.2625
## 3rd Qu.:0 3rd Qu.:-0.1297
## Max. :0 Max. : 0.9892
## NA's :2659
## tm Lm LmSex TypeLm
## Min. : 0.200 Min. : 1.50 Length:5517 Length:5517
## 1st Qu.: 2.000 1st Qu.: 16.50 Class :character Class :character
## Median : 4.000 Median : 27.10 Mode :character Mode :character
## Mean : 6.315 Mean : 41.77
## 3rd Qu.: 7.200 3rd Qu.: 49.80
## Max. :156.000 Max. :922.00
## NA's :5340 NA's :4766
## LmMale LmFemale GrowthEnviron
## Min. : 6.00 Min. : 8.00 Length:5517
## 1st Qu.: 19.85 1st Qu.: 20.40 Class :character
## Median : 29.00 Median : 34.70 Mode :character
## Mean : 48.76 Mean : 49.44
## 3rd Qu.: 50.00 3rd Qu.: 62.50
## Max. :280.00 Max. :257.00
## NA's :5410 NA's :5360There are more NAs for TLinfinity than for for K which means that for some stocks the former is missing. We thus keep only stocks with both parameters available.
nrow(growthdata[!is.na(growthdata$TLinfinity),])## [1] 5517length(unique(growthdata[!is.na(growthdata$TLinfinity),]$Species))## [1] 1313We have 5517 stocks for 1313 species with an estimate of of K and TLinfinity.
Auximetric analysis
We start by producing a log-log plot of K versus TLinfinity. According to Pauly and Munro (1984), \(log_{10}(K) \sim log_{10}(L_\infty)\) is expected to follow a regression line with a slope of 2, the intercept of which gives the growth index of the different species.
We start by adding an identifier, speciestype for 3 groups of species: tunas that are among the fish species with the higshest growth performances, deep-sea sharks that have life-history strategies relativeley close to that of coelancanths, and roughies that are among the fish with the lowest growth performances. For this, we need to load package rfishbase and use function species_list.
library(rfishbase)
growthdata$tunas<-growthdata$Species %in% species_list(Genus="Thunnus") | growthdata$Species %in% species_list(Genus="Katsuwonus") | growthdata$Species %in% species_list(Genus="Euthynnus") | growthdata$Species %in% species_list(Genus="Auxis") | growthdata$Species %in% species_list(Genus="Allothunnus")
growthdata$deepseasharks<- growthdata$Species %in% species_list(Order="Squaliformes")
growthdata$roughies<-growthdata$Species %in% species_list(Family="Trachichthyidae")
growthdata$speciestype[growthdata$tunas==T]<-"Tunas"
growthdata$speciestype[growthdata$deepseasharks==T]<-"Deep-sea sharks"
growthdata$speciestype[growthdata$roughies==T]<-"Roughies"
growthdata$speciestype[growthdata$tunas==F & growthdata$deepseasharks==F & growthdata$roughies==F]<-"Other species"
growthdata$speciestype<-as.factor(growthdata$speciestype)
summary(growthdata)## AutoCtr SpecCode Species StockCode
## Min. : 139 Min. : 4 Length:5517 Min. : 4
## 1st Qu.: 2177 1st Qu.: 204 Class :character 1st Qu.: 408
## Median : 4703 Median : 1288 Mode :character Median : 1672
## Mean : 7815 Mean : 3828 Mean : 7248
## 3rd Qu.:10235 3rd Qu.: 4449 3rd Qu.: 4986
## Max. :23914 Max. :67314 Max. :58558
##
## Locality PopGrowthRef Sex Loo
## Length:5517 Min. : 9 Length:5517 Min. : 2.10
## Class :character 1st Qu.: 2081 Class :character 1st Qu.: 25.90
## Mode :character Median : 13304 Mode :character Median : 41.00
## Mean : 38231 Mean : 62.46
## 3rd Qu.: 74354 3rd Qu.: 77.20
## Max. :119692 Max. :1966.00
##
## TLinfinity K to Winfinity
## Min. : 2.10 Min. : 0.0070 Min. :-12.210 Min. : 0
## 1st Qu.: 26.80 1st Qu.: 0.1590 1st Qu.: -1.700 1st Qu.: 193
## Median : 41.97 Median : 0.2700 Median : -0.858 Median : 721
## Mean : 63.92 Mean : 0.4517 Mean : -1.200 Mean : 31059
## 3rd Qu.: 79.00 3rd Qu.: 0.5300 3rd Qu.: -0.300 3rd Qu.: 4526
## Max. :1966.00 Max. :12.9570 Max. : 1.100 Max. :20000000
## NA's :2883 NA's :3091
## Type Auxim LinfLmax LogKLogLoo
## Length:5517 Length:5517 Min. :0 Min. :-0.8340
## Class :character Class :character 1st Qu.:0 1st Qu.:-0.4215
## Mode :character Mode :character Median :0 Median :-0.3149
## Mean :0 Mean :-0.2625
## 3rd Qu.:0 3rd Qu.:-0.1297
## Max. :0 Max. : 0.9892
## NA's :2659
## tm Lm LmSex TypeLm
## Min. : 0.200 Min. : 1.50 Length:5517 Length:5517
## 1st Qu.: 2.000 1st Qu.: 16.50 Class :character Class :character
## Median : 4.000 Median : 27.10 Mode :character Mode :character
## Mean : 6.315 Mean : 41.77
## 3rd Qu.: 7.200 3rd Qu.: 49.80
## Max. :156.000 Max. :922.00
## NA's :5340 NA's :4766
## LmMale LmFemale GrowthEnviron tunas
## Min. : 6.00 Min. : 8.00 Length:5517 Mode :logical
## 1st Qu.: 19.85 1st Qu.: 20.40 Class :character FALSE:5309
## Median : 29.00 Median : 34.70 Mode :character TRUE :208
## Mean : 48.76 Mean : 49.44
## 3rd Qu.: 50.00 3rd Qu.: 62.50
## Max. :280.00 Max. :257.00
## NA's :5410 NA's :5360
## deepseasharks roughies speciestype
## Mode :logical Mode :logical Deep-sea sharks: 58
## FALSE:5459 FALSE:5511 Other species :5245
## TRUE :58 TRUE :6 Roughies : 6
## Tunas : 208
##
##
## We then add our coelacanth data to Fishbase Von Bertalanffy growth data.
colnames(growthdata)## [1] "AutoCtr" "SpecCode" "Species" "StockCode"
## [5] "Locality" "PopGrowthRef" "Sex" "Loo"
## [9] "TLinfinity" "K" "to" "Winfinity"
## [13] "Type" "Auxim" "LinfLmax" "LogKLogLoo"
## [17] "tm" "Lm" "LmSex" "TypeLm"
## [21] "LmMale" "LmFemale" "GrowthEnviron" "tunas"
## [25] "deepseasharks" "roughies" "speciestype"coelacanth.data<-data.frame(AutoCtr=max(growthdata$AutoCtr)+1:3,SpecCode=rep(2063,3),Species=rep("Latimeria chalumnae",3),StockCode=rep(NA,3),Locality=rep("Comoro islands",3),PopGrowthRef=rep(NA,3),Sex=rep("mixed",3),Loo=c(203.8,296.1,286.4),TLinfinity=c(203.8,296.1,286.4),K=c(0.02,0.05,0.04),to=rep(NA,3),Winfinity=rep(NA,3),Type=rep("TL",3),Auxim=rep(NA,3),LinfLmax=rep(NA,3),LogKLogLoo=rep(NA,3),tm=rep(NA,3),Lm=rep(NA,3),LmSex=rep(NA,3),TypeLm=rep(NA,3),LmMale=rep(NA,3),LmFemale=rep(NA,3),GrowthEnviron=rep("open waters",3),tunas=rep(F,3),deepseasharks=rep(F,3),roughies=rep(F,3),speciestype=c("Coelacanth circuli","Coelacanth macrocirculi","Coelacanth previous study"))
head(growthdata)## AutoCtr SpecCode Species StockCode
## 1 139 4 Engraulis ringens 4
## 2 140 4 Engraulis ringens 4
## 3 141 4 Engraulis ringens 4
## 4 142 4 Engraulis ringens 4
## 5 143 4 Engraulis ringens 4
## 6 144 4 Engraulis ringens 4
## Locality PopGrowthRef Sex Loo TLinfinity K
## 1 Iquique 9 unsexed 16.9 16.9 1.60
## 2 Northern/central stock, 4-14<b0>S 9 unsexed 18.2 18.2 0.75
## 3 Northern/central stock, 4-14<b0>S 9 unsexed 18.5 18.5 0.87
## 4 Northern/central stock, 4-14<b0>S 9 unsexed 19.0 19.0 0.62
## 5 Iquique 9 unsexed 19.0 19.0 1.11
## 6 Northern/central stock, 4-14<b0>S 9 unsexed 19.3 19.3 0.75
## to Winfinity Type Auxim LinfLmax LogKLogLoo tm Lm LmSex TypeLm
## 1 -0.02 34 TL within ellipse 0 0.16624 NA NA <NA> <NA>
## 2 NA 41 TL within ellipse 0 -0.09915 NA NA <NA> <NA>
## 3 NA 43 TL within ellipse 0 -0.04773 NA NA <NA> <NA>
## 4 NA 46 TL within ellipse 0 -0.16235 NA NA <NA> <NA>
## 5 NA NA TL within ellipse 0 0.03544 NA NA <NA> <NA>
## 6 NA 48 TL within ellipse 0 -0.09719 NA NA <NA> <NA>
## LmMale LmFemale GrowthEnviron tunas deepseasharks roughies speciestype
## 1 NA NA open waters FALSE FALSE FALSE Other species
## 2 NA NA open waters FALSE FALSE FALSE Other species
## 3 NA NA open waters FALSE FALSE FALSE Other species
## 4 NA NA open waters FALSE FALSE FALSE Other species
## 5 NA NA open waters FALSE FALSE FALSE Other species
## 6 NA NA open waters FALSE FALSE FALSE Other speciesdim(growthdata)## [1] 5517 27growthdata<-rbind(growthdata,coelacanth.data)
dim(growthdata)## [1] 5520 27We can then compute the growth performance index \(\phi^\prime\) as proposed by Pauly and Munro (1984).
growthdata$gindex<-log10(growthdata$K)+2*log10(growthdata$TLinfinity)#Compute growth index log10(K)+2log10(Linfinity) for all stocks
head(growthdata)## AutoCtr SpecCode Species StockCode
## 1 139 4 Engraulis ringens 4
## 2 140 4 Engraulis ringens 4
## 3 141 4 Engraulis ringens 4
## 4 142 4 Engraulis ringens 4
## 5 143 4 Engraulis ringens 4
## 6 144 4 Engraulis ringens 4
## Locality PopGrowthRef Sex Loo TLinfinity K
## 1 Iquique 9 unsexed 16.9 16.9 1.60
## 2 Northern/central stock, 4-14<b0>S 9 unsexed 18.2 18.2 0.75
## 3 Northern/central stock, 4-14<b0>S 9 unsexed 18.5 18.5 0.87
## 4 Northern/central stock, 4-14<b0>S 9 unsexed 19.0 19.0 0.62
## 5 Iquique 9 unsexed 19.0 19.0 1.11
## 6 Northern/central stock, 4-14<b0>S 9 unsexed 19.3 19.3 0.75
## to Winfinity Type Auxim LinfLmax LogKLogLoo tm Lm LmSex TypeLm
## 1 -0.02 34 TL within ellipse 0 0.16624 NA NA <NA> <NA>
## 2 NA 41 TL within ellipse 0 -0.09915 NA NA <NA> <NA>
## 3 NA 43 TL within ellipse 0 -0.04773 NA NA <NA> <NA>
## 4 NA 46 TL within ellipse 0 -0.16235 NA NA <NA> <NA>
## 5 NA NA TL within ellipse 0 0.03544 NA NA <NA> <NA>
## 6 NA 48 TL within ellipse 0 -0.09719 NA NA <NA> <NA>
## LmMale LmFemale GrowthEnviron tunas deepseasharks roughies speciestype
## 1 NA NA open waters FALSE FALSE FALSE Other species
## 2 NA NA open waters FALSE FALSE FALSE Other species
## 3 NA NA open waters FALSE FALSE FALSE Other species
## 4 NA NA open waters FALSE FALSE FALSE Other species
## 5 NA NA open waters FALSE FALSE FALSE Other species
## 6 NA NA open waters FALSE FALSE FALSE Other species
## gindex
## 1 2.659893
## 2 2.395204
## 3 2.473863
## 4 2.349899
## 5 2.602830
## 6 2.446176We can now proceed with auximetric analysis per se. We start by producing the grid of regression lines \(log_{10}(K)\sim log_{10}(L_\infty)\) and associated labels corresponding to various intercepts \(\phi^\prime\) indicating various growth performance indices.
phitextposition<-function (phi,Klim,TLlim){
phi<-phi+0.15
Klim<-log10(Klim)
TLlim<-log10(TLlim)
(xup<-(phi-Klim[2])/2)
if (TLlim[1]<xup & xup<TLlim[2]) {
(positionLU<-10^c(xup,Klim[2]))
} else if (xup<TLlim[1]) {
yleft<-phi-2*TLlim[1]
(positionLU<-10^c(TLlim[1],yleft))
}
(xdown<-(phi-Klim[1])/2)
if (TLlim[1]<xdown & xdown<TLlim[2]) {
(positionRD<-10^c(xdown,Klim[1]))
} else if (xdown>TLlim[2]) {
yright<-phi-2*TLlim[2]
(positionRD<-10^c(TLlim[2],yright))
}
return(list(posLU=positionLU,posRD=positionRD))
}
phiindex<-0:7
Klim<-c(0.005,20)
TLlim<-c(2,2000)
gindexlabels<-lapply(phiindex,phitextposition,Klim=Klim,TLlim=TLlim)
gindexLU_xpos<-sapply(lapply(gindexlabels,'[[',1),'[[',1)
gindexLU_ypos<-sapply(lapply(gindexlabels,'[[',1),'[[',2)
gindexRD_xpos<-sapply(lapply(gindexlabels,'[[',2),'[[',1)
gindexRD_ypos<-sapply(lapply(gindexlabels,'[[',2),'[[',2)We then produce the auximetric plot using packages ggplot2 for producing plots, scales for scaling plot axes and extrafont for importing system fonts.
library(ggplot2)
library(scales)
library(extrafont)## Warning: package 'extrafont' was built under R version 4.0.3loadfonts() #loading system fonts
s<-8 #setting size of points representing coelacanth data
p1<-ggplot(data=subset(growthdata,speciestype="Other species"), aes(x = TLinfinity, y = K)) +
geom_abline(intercept=phiindex,slope=-2,linetype="dashed") + #Growth performance index Phi isolines
geom_point(aes(fill="Other species",shape="Other species"),color="white",size=6)+ #Data points for other species
geom_point(data=subset(growthdata, tunas==T),aes(fill="Tunas",shape="Tunas"),color="white",size=6)+ #Data points for Tunas
stat_ellipse(data=subset(growthdata, tunas==T),type="norm",na.rm=T,color="#CD9600",size=1.5)+ #Ellipsis for Tunas
geom_point(data=subset(growthdata, deepseasharks==T),aes(fill="Deep-sea sharks",shape="Deep-sea sharks"),color="white",size=6)+ #Data points for Deep-sea sharks
stat_ellipse(data=subset(growthdata, deepseasharks==T),type="norm",na.rm=T,color="#C77CFF",size=1.5)+ #Ellipsis for Deep-sea sharks
geom_point(data=subset(growthdata, roughies==T),aes(fill="Roughies",shape="Roughies"),color="white",size=6)+ #Data points for Roughies
stat_ellipse(data=subset(growthdata, roughies==T),type="norm",na.rm=T,color="#7CAE00",size=1.5)+ #Ellipsis for Roughies
geom_point(data=subset(growthdata, speciestype=="Coelacanth circuli"),aes(fill="Coelacanth circuli",shape="Coelacanth circuli"),size=s,color="white")+ #Data point for coelacanth using circuli for ageing; present study
geom_text(data=subset(growthdata, speciestype=="Coelacanth circuli"),aes(label="Coelacanth circuli"),hjust=-0.1,size=5)+
geom_point(data=subset(growthdata, speciestype=="Coelacanth macrocirculi"),aes(fill="Coelacanth macrocirculi",shape="Coelacanth macrocirculi"),size=s,color="white")+ #Data point for coelacanth using macro-circuli for ageing; present study
geom_text(data=subset(growthdata, speciestype=="Coelacanth macrocirculi"),aes(label="Coelacanth macrocirculi"),hjust=-0.1,size=5)+
geom_point(data=subset(growthdata, speciestype=="Coelacanth previous study"),aes(fill="Coelacanth previous study",shape="Coelacanth previous study"),size=s,color="white")+ #Data point for coelacanth using macro-circuli for ageing; previous study by Froese & Palomares
geom_text(data=subset(growthdata, speciestype=="Coelacanth previous study"),aes(label="Coelacanth previous study"),hjust=-0.1,size=5)+
geom_text(data=data.frame(phiindex,phiindex,label=paste("phi*minute == ",phiindex)),aes(label=label, x=gindexLU_xpos ,y=gindexLU_ypos),vjust=2,hjust=-0.5,parse=T,size=5)+
scale_shape_manual(name="Taxonomic\ngroup",values=c(
"Coelacanth circuli"=23,
"Coelacanth macrocirculi"=23,
"Coelacanth previous study"=23,
"Deep-sea sharks"=21,
"Other species"=21,
"Roughies"=21,
"Tunas"=21))+
scale_fill_manual(name="Taxonomic\ngroup",values=c(
"Coelacanth circuli"="#F8766D",
"Coelacanth macrocirculi"="#00BFC4",
"Coelacanth previous study"="#619CFF",
"Deep-sea sharks"="#C77CFF",
"Other species"="darkgrey",
"Roughies"="#7CAE00",
"Tunas"="#CD9600"))+
theme(legend.position="left",
legend.key=element_blank(),
legend.title=element_text(size=14),
legend.text=element_text(size=12),
axis.title.x = element_text(size=16, face="bold"),
axis.title.y = element_text(size=16, face="bold"),
axis.text=element_text(size=12),
text=element_text(family="Helvetica"))+
xlab(expression(bold("Asymptotic total length,"~TL[infinity]~(cm)))) +
ylab(expression(bold("Growth rate coefficient,"~K~(year^-1))))+
scale_x_log10(breaks = trans_breaks("log10", function(x) 10^x,n=3),
labels = trans_format("log10", math_format(10^.x)),limits=c(2,2000),expand=c(0,0))+
scale_y_log10(breaks = trans_breaks("log10", function(x) 10^x),
labels = trans_format("log10", math_format(10^.x)),limits=c(0.005,20),expand=c(0,0))+
annotation_logticks(sides="trbl")+
guides(shape = guide_legend(override.aes = list(size = 5)))In addition to the auximetric plot, we produce a simple boxplot of the growth performance index \(\phi^\prime\) of the different species types (tunas, deep-sea sharks, roughies) to position coelacanth estimates along this univariate dimension.
We start by computing \(\phi^\prime\) for the different coelacanth estimates
gindex_circuli<-growthdata[growthdata$speciestype=="Coelacanth circuli",]$gindex
gindex_macrocirculi_pst<-growthdata[growthdata$speciestype=="Coelacanth macrocirculi",]$gindex
gindex_macrocirculi_FroesePalomares<-growthdata[growthdata$speciestype=="Coelacanth previous study",]$gindex
gindex_coelacanth<-c(gindex_circuli,gindex_macrocirculi_pst,gindex_macrocirculi_FroesePalomares)
quantilelabels<-sapply(growthdata[growthdata$speciestype %in% c("Coelacanth circuli","Coelacanth macrocirculi","Coelacanth previous study"),]$gindex,function (g) round(100*ecdf(growthdata$gindex)(g),0))Then we produce the boxplot
p2<-ggplot(data=subset(growthdata,!(speciestype %in% c("Coelacanth circuli","Coelacanth macrocirculi","Coelacanth previous study"))), aes(x=speciestype, y=gindex)) +
geom_boxplot(aes(fill=speciestype),show.legend=F,outlier.size=3)+ #Boxplot per species type
stat_summary(fun=mean, geom="point", shape=21, fill="black",size=4) +
stat_summary(data=data.frame(speciestype="Other species",gindex=gindex_circuli),fun=mean, geom="point", shape=23,fill="#F8766D",color="white", size=6) + #Adding data poitn for coelacanth aged using circuli; present study
geom_hline(aes(yintercept=gindex_circuli),color="#F8766D") +
stat_summary(data=data.frame(speciestype="Other species",gindex=gindex_macrocirculi_pst),fun=mean, geom="point", shape=23,fill="#00BFC4",color="white", size=6) + #Adding data poitn for coelacanth aged using macro-circuli; present study
geom_hline(aes(yintercept=gindex_macrocirculi_pst),color="#00BFC4") +
stat_summary(data=data.frame(speciestype="Other species",gindex=gindex_macrocirculi_FroesePalomares),fun=mean, geom="point", shape=23,fill="#619CFF",color="white", size=6) + #Adding data poitn for coelacanth aged using macro-circuli; previous study by Froese & Palomares
geom_hline(aes(yintercept=gindex_macrocirculi_FroesePalomares),color="#619CFF") +
geom_text(data=data.frame(label=paste(quantilelabels,"%",sep="")),aes(label=label, x=-Inf ,y=gindex_coelacanth),vjust=-0.25,size=5)+
scale_fill_manual(values=c(
"Deep-sea sharks"="#C77CFF",
"Other species"="darkgrey",
"Roughies"="#7CAE00",
"Tunas"="#CD9600"))+
labs(x="Taxonomic group", y=expression(bold("Growth index,"~phi*minute==log[10](K)+2*Log[10](TL[infinity])))) +
coord_flip()+
theme(
axis.title.x = element_text(size=16, face="bold"),
axis.title.y = element_text(size=16, face="bold"),
axis.text=element_text(size=12),
legend.key = element_blank(),
text=element_text(family="Helvetica")#, legend.box = "horizontal"
)We finally combine the two plots using package ggpubr and save them.
library(ggpubr)
p_auxi<-ggarrange(p1,p2,labels=c("A","B"),ncol=1,nrow=2,heights=c(2,1),common.legend=T,legend="left")
p_auxi
ggsave("Auximetric plot_Helvetica.pdf",p_auxi,width=45,height=27,units="cm")Slow-fast life-history continuum analysis
Life-history data
We start importing the data. These data were extracted from Fishbase using the package rfishbase on June 2020. Only species from brackish and salted waters were included.
Von Bertalanffy growth parameters were extracted as described in the previous section.
Maturity data, i.e. average age at maturity \(t_m\) and length at maturity \(L_m\), were extracted from the Maturity table. For stocks for whicha range of ages and/or lengths at maturity were available instead of a point estimate, we computed the mid-age and/or -length and used them as estimates of \(t_m\) and/or \(L_m\). If not a range but a minimum or maximum age or length at maturity was provided, we used it directly as an estimate for \(t_m\) or \(L_m\). We then extracted stocks for which we had at least an estimate of either \(t_m\) or \(L_m\), as \(t_m\) could be estimated from \(L_m\) by inverting the Von Bertalanffy growth model when available. When \(L_m\) was not given as total length (or fork length for Scombridae), it was converted using length-length relationships from the Length-Length table when the information was available.
Gestation time or spawning frequency data were extracted from the Fecundity table cross checked with the Spawning table.
Maximum longevity was extracted from the Species table.
Data were merged in a single table containing life-history traits with time dimensionality, namely Von Bertalanffy rate coefficient K, age at maturity tm, maximum longevity LongevityWild and spawning frequency SpawningFreq, plus asymptotic total length TLinfinity for correcting for size effects (see below).
LHdata<-read.csv("Life_history_traits_Fishbase_2020_June.csv")
summary(LHdata)## Species SpecCode StockCode Sex
## Length:4443 Min. : 6.0 Min. : 38 Length:4443
## Class :character 1st Qu.: 30.0 1st Qu.: 85 Class :character
## Mode :character Median : 142.0 Median : 520 Mode :character
## Mean : 945.3 Mean : 7445
## 3rd Qu.: 1327.0 3rd Qu.: 1948
## Max. :16324.0 Max. :56833
##
## TLinfinity K to tunas
## Min. : 2.10 Min. : 0.0200 Min. :-8.7800 Mode :logical
## 1st Qu.: 32.60 1st Qu.: 0.1390 1st Qu.:-1.0120 FALSE:3946
## Median : 54.70 Median : 0.2400 Median :-0.0300 TRUE :497
## Mean : 77.57 Mean : 0.3475 Mean :-0.6913
## 3rd Qu.: 88.70 3rd Qu.: 0.4400 3rd Qu.: 0.0000
## Max. :651.00 Max. :12.9570 Max. : 0.7850
##
## deepseasharks roughies speciestype gindex
## Mode :logical Mode :logical Length:4443 Min. :1.402
## FALSE:4435 FALSE:4430 Class :character 1st Qu.:2.547
## TRUE :8 TRUE :13 Mode :character Median :2.826
## Mean :2.916
## 3rd Qu.:3.225
## Max. :4.788
##
## TLm tm LongevityWild SpawningFreq
## Min. : 1.0 Min. : 0.082 Min. : 0.16 Min. : 0.330
## 1st Qu.: 20.2 1st Qu.: 1.500 1st Qu.: 9.00 1st Qu.: 1.000
## Median : 30.6 Median : 3.000 Median : 20.00 Median : 1.000
## Mean : 45.8 Mean : 4.156 Mean : 20.94 Mean : 1.856
## 3rd Qu.: 43.5 3rd Qu.: 5.012 3rd Qu.: 25.00 3rd Qu.: 1.000
## Max. :427.4 Max. :86.949 Max. :149.00 Max. :150.000
## NA's :1432Then we add coelacanth data.
coelacanth.data[,c("Species","SpecCode","StockCode","Sex","TLinfinity","K","to","tm","deepseasharks","roughies","speciestype")]## Species SpecCode StockCode Sex TLinfinity K to tm
## 1 Latimeria chalumnae 2063 NA mixed 203.8 0.02 NA NA
## 2 Latimeria chalumnae 2063 NA mixed 296.1 0.05 NA NA
## 3 Latimeria chalumnae 2063 NA mixed 286.4 0.04 NA NA
## deepseasharks roughies speciestype
## 1 FALSE FALSE Coelacanth circuli
## 2 FALSE FALSE Coelacanth macrocirculi
## 3 FALSE FALSE Coelacanth previous study(combdata.ceolacanth<-growthdata[growthdata$Species=="Latimeria chalumnae",c("Species","SpecCode","StockCode","Sex","TLinfinity","K","to","tunas","deepseasharks","roughies","speciestype","gindex")])## Species SpecCode StockCode Sex TLinfinity K to tunas
## 3741 Latimeria chalumnae 2063 2258 mixed 219.0 0.058 -3.6 FALSE
## 5518 Latimeria chalumnae 2063 NA mixed 203.8 0.020 NA FALSE
## 5519 Latimeria chalumnae 2063 NA mixed 296.1 0.050 NA FALSE
## 5520 Latimeria chalumnae 2063 NA mixed 286.4 0.040 NA FALSE
## deepseasharks roughies speciestype gindex
## 3741 FALSE FALSE Other species 3.444316
## 5518 FALSE FALSE Coelacanth circuli 2.919438
## 5519 FALSE FALSE Coelacanth macrocirculi 3.641847
## 5520 FALSE FALSE Coelacanth previous study 3.516006(combdata.ceolacanth<-combdata.ceolacanth[-1,])## Species SpecCode StockCode Sex TLinfinity K to tunas
## 5518 Latimeria chalumnae 2063 NA mixed 203.8 0.02 NA FALSE
## 5519 Latimeria chalumnae 2063 NA mixed 296.1 0.05 NA FALSE
## 5520 Latimeria chalumnae 2063 NA mixed 286.4 0.04 NA FALSE
## deepseasharks roughies speciestype gindex
## 5518 FALSE FALSE Coelacanth circuli 2.919438
## 5519 FALSE FALSE Coelacanth macrocirculi 3.641847
## 5520 FALSE FALSE Coelacanth previous study 3.516006combdata.ceolacanth$TLm<-147 #Median maturation length of males and females
combdata.ceolacanth$TL1<-c(20.14,46.91,21.98)
combdata.ceolacanth$tm<-with(combdata.ceolacanth,1-(1/K)*log((TLinfinity-TLm)/(TLinfinity-TL1)))
combdata.ceolacanth$LongevityWild<-c(85,17,20)
combdata.ceolacanth$SpawningFreq<-c(1/5,1/1.5,1/1.5)
combdata.ceolacanth<-combdata.ceolacanth[,-which(colnames(combdata.ceolacanth)=="TL1")]
head(LHdata)## Species SpecCode StockCode Sex TLinfinity K to tunas
## 1 Achoerodus gouldii 14420 14219 female 68.2 0.140 0.06 FALSE
## 2 Alopias pelagicus 5891 6200 mixed 328.0 0.120 0.00 FALSE
## 3 Alopias superciliosus 2534 2730 male 385.0 0.088 -4.24 FALSE
## 4 Alopias superciliosus 2534 2730 female 422.0 0.092 -4.21 FALSE
## 5 Alopias superciliosus 2534 2730 female 422.0 0.092 -4.21 FALSE
## 6 Alopias superciliosus 2534 2730 male 385.0 0.088 -4.24 FALSE
## deepseasharks roughies speciestype gindex TLm tm LongevityWild
## 1 FALSE FALSE Other species 2.813697 65.30 17.00000 70
## 2 FALSE FALSE Other species 4.110929 NA 13.20000 29
## 3 FALSE FALSE Other species 4.115404 276.00 18.57972 20
## 4 FALSE FALSE Other species 4.214413 NA 11.50000 20
## 5 FALSE FALSE Other species 4.214413 374.96 28.05790 20
## 6 FALSE FALSE Other species 4.115404 NA 10.00000 20
## SpawningFreq
## 1 1
## 2 1
## 3 1
## 4 1
## 5 1
## 6 1dim(LHdata)## [1] 4443 16LHdata<-rbind(LHdata,combdata.ceolacanth)
dim(LHdata)## [1] 4446 16Slow-fast continuum analysis
Data inspection
We start by looking at the data after \(log_{10}\) transformation
panel.hist <- function(x, ...)
{
usr <- par("usr"); on.exit(par(usr))
par(usr = c(usr[1:2], 0, 1.5) )
h <- hist(x, plot = FALSE)
breaks <- h$breaks; nB <- length(breaks)
y <- h$counts; y <- y/max(y)
rect(breaks[-nB], 0, breaks[-1], y, col = "cyan", ...)
}
panel.cor <- function(x, y, digits = 2, prefix = "", cex.cor, ...)
{
usr <- par("usr"); on.exit(par(usr))
par(usr = c(0, 1, 0, 1))
r <- abs(cor(x, y,use="complete.obs"))
txt <- format(c(r, 0.123456789), digits = digits)[1]
txt <- paste0(prefix, txt)
if(missing(cex.cor)) cex.cor <- 0.8/strwidth(txt)
text(0.5, 0.5, txt, cex = cex.cor * r)
}
pairs(log10(LHdata[c("TLinfinity","K","gindex","tm","LongevityWild","SpawningFreq")]),diag.panel = panel.hist,lower.panel=panel.smooth,upper.panel=panel.cor)
head(LHdata)## Species SpecCode StockCode Sex TLinfinity K to tunas
## 1 Achoerodus gouldii 14420 14219 female 68.2 0.140 0.06 FALSE
## 2 Alopias pelagicus 5891 6200 mixed 328.0 0.120 0.00 FALSE
## 3 Alopias superciliosus 2534 2730 male 385.0 0.088 -4.24 FALSE
## 4 Alopias superciliosus 2534 2730 female 422.0 0.092 -4.21 FALSE
## 5 Alopias superciliosus 2534 2730 female 422.0 0.092 -4.21 FALSE
## 6 Alopias superciliosus 2534 2730 male 385.0 0.088 -4.24 FALSE
## deepseasharks roughies speciestype gindex TLm tm LongevityWild
## 1 FALSE FALSE Other species 2.813697 65.30 17.00000 70
## 2 FALSE FALSE Other species 4.110929 NA 13.20000 29
## 3 FALSE FALSE Other species 4.115404 276.00 18.57972 20
## 4 FALSE FALSE Other species 4.214413 NA 11.50000 20
## 5 FALSE FALSE Other species 4.214413 374.96 28.05790 20
## 6 FALSE FALSE Other species 4.115404 NA 10.00000 20
## SpawningFreq
## 1 1
## 2 1
## 3 1
## 4 1
## 5 1
## 6 1We can already see that the growth index \(\phi^\prime\) is very much determined by \(TL_\infty\). There is also a strong correlation between the von Bertalanffy rate \(K\), age at maturity \(t_m\) and longevity which is quite expected in case of a slow-fast continuum. The correlation with spawning cycles is less obvious.
Analysis without accounting for phylogeny
We proceeed with the first analysis that looks at the slow-fast continuum after removing the size effect but wihtout accounting for phylogenetic relationships bewteen species. For this, we regress each log10-transformed life-history traits with time dimensionality on log10-transformed \(TL_\infty\) and take the residuals. We then apply a centred and scaled PCA on the residuals. We use packages ade4 and factoextra.
library(ade4)
library(factoextra) #load package for ggplot-like plotting of pca results## Welcome! Want to learn more? See two factoextra-related books at https://goo.gl/ve3WBapartial.matrix<-sapply(LHdata[c("K","tm","LongevityWild","SpawningFreq")],function (x) resid(lm(log10(x)~log10(LHdata$TLinfinity)))) #Partialling out size effect through linear regression
pc.cor.partial<-dudi.pca(partial.matrix, scannf = F,nf=5) #Perfoming centred and scaled PCA
(var.expl.partial<-round(100*pc.cor.partial$eig/sum(pc.cor.partial$eig),1))## [1] 54.6 24.5 12.5 8.5cumsum(var.expl.partial) #The 2 first axes explain 79.1%## [1] 54.6 79.1 91.6 100.1We produce the customized biplot of the PCA.
pcresults.partial<-data.frame(pc.cor.partial$li,speciestype=LHdata$speciestype) #Gathering results
s<-8 #setting size of points representing coelacanth data
p3<-ggplot(data=subset(pcresults.partial,speciestype=="Other species"),aes(x=Axis1,y=Axis2))+
geom_point(aes(fill="Other species",shape="Other species"),color="white",size=6)+ #Other species stocks projection on the plane of the 2 first PCs
geom_point(data=subset(pcresults.partial,speciestype=="Tunas"),aes(fill="Tunas",shape="Tunas"),color="white",size=6)+ #Tunas' stocks projection on the plane of the 2 first PCs
stat_ellipse(data=subset(pcresults.partial,speciestype=="Tunas"),type="norm",na.rm=T,color="#CD9600",size=1)+ #Ellipsis for tunas
geom_point(data=subset(pcresults.partial,speciestype=="Deep-sea sharks"),aes(fill="Deep-sea sharks",shape="Deep-sea sharks"),color="white",size=6)+ #Deep-sea sharks' stocks projection on the plane of the 2 first PCs
stat_ellipse(data=subset(pcresults.partial,speciestype=="Deep-sea sharks"),type="norm",na.rm=T,color="#C77CFF",size=1)+ #Ellipsis for deep-sea sharks
geom_point(data=subset(pcresults.partial,speciestype=="Roughies"),aes(fill="Roughies",shape="Roughies"),color="white",size=6)+ #Roughies' stocks projection on the plane of the 2 first PCs
stat_ellipse(data=subset(pcresults.partial,speciestype=="Roughies"),type="norm",na.rm=T,color="#7CAE00",size=1)+ #Ellipsis for roughies
geom_point(data=subset(pcresults.partial,speciestype=="Coelacanth circuli"),aes(fill="Coelacanth circuli",shape="Coelacanth circuli"),size=s,color="white")+ #Coelacanth projection when aged using circuli; present study
geom_text(data=subset(pcresults.partial, speciestype=="Coelacanth circuli"),aes(label="Coelacanth circuli"),hjust=0.5,vjust=-1.4,size=5)+
geom_point(data=subset(pcresults.partial,speciestype=="Coelacanth macrocirculi"),aes(fill="Coelacanth macrocirculi",shape="Coelacanth macrocirculi"),size=s,color="white")+ #Coelacanth projection when aged using macro-circuli; present study
geom_text(data=subset(pcresults.partial, speciestype=="Coelacanth macrocirculi"),aes(label="Coelacanth macrocirculi"),hjust=0.1,vjust=2,size=5)+
geom_point(data=subset(pcresults.partial,speciestype=="Coelacanth previous study"),aes(fill="Coelacanth previous study",shape="Coelacanth previous study"),size=s,color="white")+ #Coelacanth projection when aged using macro-circuli; previosu study by Froese & Palomares
geom_text(data=subset(pcresults.partial, speciestype=="Coelacanth previous study"),aes(label="Coelacanth previous study"),hjust=0.85,vjust=2,size=5)+
geom_segment(data=pc.cor.partial$c1,aes(x=rep(0,4),y=rep(0,4),xend=5*CS1,yend=5*CS2),arrow=arrow(length=unit(0.015, "npc")),size=1)+ #Arrows representing life-history traits
geom_text(data=pc.cor.partial$c1,aes(x=5*CS1,y=5*CS2,label=c("Growth coefficient","Maturation age","Longevity","Spawning frequency")),nudge_x=0.75*c(1.2,-0.9,-0.6,0),nudge_y=0.2*c(0,-0.15,1,1),size=5)+
scale_shape_manual(name="Taxonomic\ngroup",values=c(
"Coelacanth circuli"=23,
"Coelacanth macrocirculi"=23,
"Coelacanth previous study"=23,
"Deep-sea sharks"=21,
"Other species"=21,
"Roughies"=21,
"Tunas"=21))+
scale_fill_manual(name="Taxonomic\ngroup",values=c(
"Coelacanth circuli"="#F8766D",
"Coelacanth macrocirculi"="#00BFC4",
"Coelacanth previous study"="#619CFF",
"Deep-sea sharks"="#C77CFF",
"Other species"="darkgrey",
"Roughies"="#7CAE00",
"Tunas"="#CD9600"))+
xlab(paste("PC1 (",var.expl.partial[1],"%)",sep=""))+
ylab(paste("PC2 (",var.expl.partial[2],"%)",sep=""))+
scale_x_continuous(breaks=c(-4,-2,0,2,4))+
scale_y_continuous(breaks=c(-2.5,0,2.5,5,7.5))+
theme(legend.position="left",
legend.key=element_blank(),
legend.title=element_text(size=14),
legend.text=element_text(size=12),
axis.title.x = element_text(size=16, face="bold"),
axis.title.y = element_text(size=16, face="bold"),
axis.text=element_text(size=12),
text=element_text(family="Helvetica"))+
geom_hline(yintercept=0,linetype=2)+
geom_vline(xintercept=0,linetype=2)+
guides(shape = guide_legend(override.aes = list(size = 5)))In addition to the PCA biplot, we produce a simple boxplot of the stock scores along the first PC taken an index of the slow-fast continuum for the different species types (tunas, deep-sea sharks, roughies) to position coelacanth estimates along this univariate dimension.
We start by extracting the scores for the different coelacanth estimates
Axis1_circuli<-pcresults.partial[pcresults.partial$speciestype=="Coelacanth circuli",]$Axis1
Axis1_macrocirculi_pst<-pcresults.partial[pcresults.partial$speciestype=="Coelacanth macrocirculi",]$Axis1
Axis1_macrocirculi_FroesePalomares<-pcresults.partial[pcresults.partial$speciestype=="Coelacanth previous study",]$Axis1
Axis1_coelacanth<-c(Axis1_circuli,Axis1_macrocirculi_pst,Axis1_macrocirculi_FroesePalomares)
quantilelabels.sf<-sapply(pcresults.partial[pcresults.partial$speciestype %in% c("Coelacanth circuli","Coelacanth macrocirculi","Coelacanth previous study"),]$Axis1,function (g) round(100*ecdf(pcresults.partial$Axis1)(g),0))Then, we produce the boxplot.
p4<-ggplot(data=subset(pcresults.partial,!(speciestype %in% c("Coelacanth circuli","Coelacanth macrocirculi","Coelacanth previous study"))), aes(x=speciestype, y=Axis1)) +
geom_boxplot(aes(fill=speciestype),show.legend=F,outlier.size=3)+ #Boxplot per species type
stat_summary(fun=mean, geom="point", shape=21, fill="black",size=4) +
stat_summary(data=data.frame(speciestype="Other species",Axis1=Axis1_circuli),fun=mean, geom="point", shape=23,fill="#F8766D",color="white", size=6) + #Adding coelacanth aged using circuli; present study
geom_hline(aes(yintercept=Axis1_circuli),color="#C77CFF") +
stat_summary(data=data.frame(speciestype="Other species",Axis1=Axis1_macrocirculi_pst),fun=mean, geom="point", shape=23,fill="#00BFC4",color="white", size=6) + #Adding coelacanth aged using macro-circuli; present study
geom_hline(aes(yintercept=Axis1_macrocirculi_pst),color="#CD9600") +
stat_summary(data=data.frame(speciestype="Other species",Axis1=Axis1_macrocirculi_FroesePalomares),fun=mean, geom="point", shape=23,fill="#619CFF",color="white", size=6) + #Adding coelacanth aged using macro-circuli; previous study by Froese & Palomares
geom_hline(aes(yintercept=Axis1_macrocirculi_FroesePalomares),color="#619CFF") +
geom_text(data=data.frame(label=paste(quantilelabels.sf,"%",sep="")),aes(label=label, x=-Inf ,y=Axis1_coelacanth),vjust=-0.25,size=5)+
scale_fill_manual(values=c(
"Deep-sea sharks"="#C77CFF",
"Other species"="darkgrey",
"Roughies"="#7CAE00",
"Tunas"="#CD9600"))+
labs(x="Taxonomic group", y="Slow-fast life-history axis, PC1") +
coord_flip()+
theme(
axis.title.x = element_text(size=16, face="bold"),
axis.title.y = element_text(size=16, face="bold"),
axis.text=element_text(size=12),
legend.key = element_blank(),
text=element_text(family="Helvetica")
)We finally combine the two plots and save them.
p_pca<-ggarrange(p3,p4,labels=c("A","B"),ncol=1,nrow=2,heights=c(2,1),common.legend=T,legend="left")
p_pca
ggsave("Slow-fast life-history_Helvetica.pdf",p_pca,width=45,height=27,units="cm")Analysis accounting for phylogeny
Then we analyze the life-history data accounting for phylogenetic relationships between species. The problem is no fine-grained (i.e. including all species) phylogeney relating teleosts,elasmobranchs and coealcanth is available. We perfomr the analysis only on teleost as a sufciently fine-grained phylogney is available for this group to check whether the slow-fast continuum is modified. For this we use the Fish Tree of Life phylogeny that can be extrated using the package fishtree.
We strat extracting the phylogeny of species in the life-history dataset.
summary(LHdata)## Species SpecCode StockCode Sex
## Length:4446 Min. : 6 Min. : 38 Length:4446
## Class :character 1st Qu.: 30 1st Qu.: 85 Class :character
## Mode :character Median : 142 Median : 520 Mode :character
## Mean : 946 Mean : 7445
## 3rd Qu.: 1327 3rd Qu.: 1948
## Max. :16324 Max. :56833
## NA's :3
## TLinfinity K to tunas
## Min. : 2.10 Min. : 0.0200 Min. :-8.7800 Mode :logical
## 1st Qu.: 32.60 1st Qu.: 0.1390 1st Qu.:-1.0120 FALSE:3949
## Median : 54.70 Median : 0.2400 Median :-0.0300 TRUE :497
## Mean : 77.69 Mean : 0.3473 Mean :-0.6913
## 3rd Qu.: 88.70 3rd Qu.: 0.4400 3rd Qu.: 0.0000
## Max. :651.00 Max. :12.9570 Max. : 0.7850
## NA's :3
## deepseasharks roughies speciestype gindex
## Mode :logical Mode :logical Length:4446 Min. :1.402
## FALSE:4438 FALSE:4433 Class :character 1st Qu.:2.547
## TRUE :8 TRUE :13 Mode :character Median :2.826
## Mean :2.916
## 3rd Qu.:3.229
## Max. :4.788
##
## TLm tm LongevityWild SpawningFreq
## Min. : 1.0 Min. : 0.082 Min. : 0.16 Min. : 0.200
## 1st Qu.: 20.2 1st Qu.: 1.500 1st Qu.: 9.00 1st Qu.: 1.000
## Median : 30.6 Median : 3.000 Median : 20.00 Median : 1.000
## Mean : 45.9 Mean : 4.173 Mean : 20.96 Mean : 1.855
## 3rd Qu.: 43.5 3rd Qu.: 5.046 3rd Qu.: 25.00 3rd Qu.: 1.000
## Max. :427.4 Max. :86.949 Max. :149.00 Max. :150.000
## NA's :1432library(fishtree)## Warning: package 'fishtree' was built under R version 4.0.5taxo<-fishtree_taxonomy()
taxo[taxo$name=="Squalus",]## [1] rank name
## <0 rows> (or 0-length row.names)str(taxo)## 'data.frame': 685 obs. of 2 variables:
## $ rank: Factor w/ 18 levels "class","subclass",..: 1 1 2 2 3 3 4 4 5 5 ...
## $ name: chr "Actinopteri" "Cladistia" "Chondrostei" "Neopterygii" ...phy <- fishtree_phylogeny(species = unique((LHdata$Species)),type='phylogram') #Importing phylogny of species present in the life-history dataset## Warning: Requested 132 but only found 105 species.
## * Achoerodus gouldii
## * Alopias pelagicus
## * Alopias superciliosus
## * Alopias vulpinus
## * Ammodytes marinus
## * ...and 22 othersplot #Plot phylo tree## function (x, y, ...)
## UseMethod("plot")
## <bytecode: 0x0000000038c323f0>
## <environment: namespace:base>(abs.species<-sub(" ","_",unique(LHdata$Species))[!(sub(" ","_",unique(LHdata$Species)) %in% phy$tip.label)]) #Identifying missing species in the phylogeny compared to life-history dataset## [1] "Achoerodus_gouldii" "Alopias_pelagicus"
## [3] "Alopias_superciliosus" "Alopias_vulpinus"
## [5] "Ammodytes_marinus" "Carcharhinus_acronotus"
## [7] "Carcharhinus_brachyurus" "Carcharhinus_falciformis"
## [9] "Carcharhinus_limbatus" "Carcharhinus_obscurus"
## [11] "Carcharhinus_plumbeus" "Carcharhinus_sorrah"
## [13] "Eviota_sigillata" "Isogomphodon_oxyrhynchus"
## [15] "Lamna_nasus" "Latimeria_chalumnae"
## [17] "Mustelus_mustelus" "Nemipterus_randalli"
## [19] "Prionace_glauca" "Rhizoprionodon_taylori"
## [21] "Rhizoprionodon_terraenovae" "Scyliorhinus_canicula"
## [23] "Sebastes_goodei" "Sphyrna_lewini"
## [25] "Sphyrna_tiburo" "Squalus_megalops"
## [27] "Triakis_semifasciata"LHdata.phylo<-LHdata[!(sub(" ","_",LHdata$Species) %in% abs.species),]#Remove Life-history data from species that are not in the phylogeny treeThen we perform a phylogenetic PCA using package phytools after correcting for size effects using phylogenetic regression.
We start by size effect correction using function pgls.Ives and group the residual sin a common data frame resid.LHdata.phylo. Total length and klife history traits with time dimensionality are log10-transformed.
library(phytools)## Warning: package 'phytools' was built under R version 4.0.5## Loading required package: ape## Warning: package 'ape' was built under R version 4.0.5##
## Attaching package: 'ape'## The following object is masked from 'package:ggpubr':
##
## rotate## Loading required package: maps## Warning: package 'maps' was built under R version 4.0.3#phyl.resid(tree=phy,x=LHdata$TLinfinity,Y=LHdata[c("K","tm","LongevityWild","SpawningFreq")])
named.TLinfinity<-log10(LHdata.phylo$TLinfinity)
names(named.TLinfinity)<-sub(" ","_",LHdata.phylo$Species)
named.K<-log10(LHdata.phylo$K) #Size-correction for K
names(named.K)<-sub(" ","_",LHdata.phylo$Species)
(K.phylo<-pgls.Ives(tree=phy,X=named.TLinfinity,y=named.K))##
## ---------------------------------------------------------------
## | **Warning: |
## | User reports suggest that this method may frequently |
## | fail to find the ML solution. Please use with caution. |
## ---------------------------------------------------------------## Warning: Some species contain only one sample. Substituting mean variance.## Warning: Some species contain only one sample. Substituting mean variance.
## Warning: Some species contain only one sample. Substituting mean variance.##
## Result PGLS with sampling error in x & y
## (based on Ives et al. 2007):
##
## Response: y
## Intercept beta[1]
## y 0.5109033 -0.7168266
##
## Summary of ML estimated parameters:
## sigma^2 a log(L)
## y 0.330662 -0.703630 -30.332793
## x 0.348500 1.694319
## ---------
##
## R thinks it has converged.resid.K<-named.K-(K.phylo$beta[1]+K.phylo$beta[2]*named.TLinfinity)
named.tm<-log10(LHdata.phylo$tm) #Size-correction for tm
names(named.tm)<-sub(" ","_",LHdata.phylo$Species)
(tm.phylo<-pgls.Ives(tree=phy,X=named.TLinfinity,y=named.tm))##
## ---------------------------------------------------------------
## | **Warning: |
## | User reports suggest that this method may frequently |
## | fail to find the ML solution. Please use with caution. |
## ---------------------------------------------------------------## Warning: Some species contain only one sample. Substituting mean variance.
## Warning: Some species contain only one sample. Substituting mean variance.
## Warning: Some species contain only one sample. Substituting mean variance.##
## Result PGLS with sampling error in x & y
## (based on Ives et al. 2007):
##
## Response: y
## Intercept beta[1]
## y -0.2662341 0.5066253
##
## Summary of ML estimated parameters:
## sigma^2 a log(L)
## y 0.390831 0.591775 -54.954259
## x 0.356336 1.693576
## ---------
##
## R thinks it has converged.resid.tm<-named.tm-(tm.phylo$beta[1]+tm.phylo$beta[2]*named.TLinfinity)
named.LongevityWild<-log10(LHdata.phylo$LongevityWild) #Size-correction for Longevity
names(named.LongevityWild)<-sub(" ","_",LHdata.phylo$Species)
(LongevityWild.phylo<-pgls.Ives(tree=phy,X=named.TLinfinity,y=named.LongevityWild))##
## ---------------------------------------------------------------
## | **Warning: |
## | User reports suggest that this method may frequently |
## | fail to find the ML solution. Please use with caution. |
## ---------------------------------------------------------------## Warning: Some species contain only one sample. Substituting mean variance.
## Warning: Some species contain only one sample. Substituting mean variance.
## Warning: Some species contain only one sample. Substituting mean variance.##
## Result PGLS with sampling error in x & y
## (based on Ives et al. 2007):
##
## Response: y
## Intercept beta[1]
## y -0.05247589 0.7365154
##
## Summary of ML estimated parameters:
## sigma^2 a log(L)
## y 0.584115 1.195408 -48.985933
## x 0.370023 1.694307
## ---------
##
## R thinks it has converged.resid.LongevityWild<-named.LongevityWild-(LongevityWild.phylo$beta[1]+LongevityWild.phylo$beta[2]*named.TLinfinity)
named.SpawningFreq<-log10(LHdata.phylo$SpawningFreq) #Size-correction for Spawning frequency
names(named.SpawningFreq)<-sub(" ","_",LHdata.phylo$Species)
(SpawningFreq.phylo<-pgls.Ives(tree=phy,X=named.TLinfinity,y=named.SpawningFreq))##
## ---------------------------------------------------------------
## | **Warning: |
## | User reports suggest that this method may frequently |
## | fail to find the ML solution. Please use with caution. |
## ---------------------------------------------------------------## Warning: Some species contain only one sample. Substituting mean variance.
## Warning: Some species contain only one sample. Substituting mean variance.
## Warning: Some species contain only one sample. Substituting mean variance.##
## Result PGLS with sampling error in x & y
## (based on Ives et al. 2007):
##
## Response: y
## Intercept beta[1]
## y 0.1760628 -0.07051426
##
## Summary of ML estimated parameters:
## sigma^2 a log(L)
## y 0.827221 0.056748 -65.124657
## x 0.349658 1.692068
## ---------
##
## R thinks it has converged.resid.SpawningFreq<-named.SpawningFreq-(SpawningFreq.phylo$beta[1]+SpawningFreq.phylo$beta[2]*named.TLinfinity)
#Gathering residuals in a single data frame
resid.LHdata.phylo<-data.frame(K=resid.K,tm=resid.tm,LongevityWild=resid.LongevityWild,SpawningFreq=resid.SpawningFreq)
resid.LHdata.phylo$Species<-LHdata.phylo$Species
resid.LHdata.phylo$speciestype<-LHdata.phylo$speciestypeThen we perform the phylogenetic PCA on the averages per species and project the single stock scores on the resulting PCA plane.
(average.resid.LHdata.phylo<-aggregate(cbind(K,tm,LongevityWild,SpawningFreq)~Species,data=resid.LHdata.phylo,FUN=mean)) #Averaging by species## Species K tm LongevityWild
## 1 Anoplopoma fimbria 0.133605686 0.036514210 0.628241328
## 2 Atractoscion nobilis 0.147772321 -0.228221897 -0.240573262
## 3 Aulorhynchus flavidus 0.136152464 -0.391797789 0.050093143
## 4 Beryx splendens -0.121965077 0.164921460 0.129075931
## 5 Caranx crysos 0.228592182 -0.012072794 -0.129275401
## 6 Centropristis striata 0.250381545 -0.277177490 0.214774788
## 7 Chelidonichthys lucerna -0.165935155 -0.179376822 -0.003329178
## 8 Clupea harengus 0.060498712 0.038451170 0.360566184
## 9 Coryphaena hippurus 1.151371387 -1.272453628 -1.002431845
## 10 Coryphaenoides rupestris -0.120869106 0.276869155 0.361465195
## 11 Cynoscion nebulosus 0.114845140 -0.286421744 -0.046178159
## 12 Diplodus annularis -0.201040215 0.033983409 -0.067629115
## 13 Engraulis encrasicolus 0.091490260 -0.251644230 -0.167626338
## 14 Engraulis mordax -0.004639927 -0.105730913 -0.045528672
## 15 Eopsetta jordani -0.089983477 0.229413680 0.332830233
## 16 Epinephelus areolatus 0.098977047 -0.074271649 0.076851961
## 17 Epinephelus chlorostigma 0.059012944 0.011153074 0.259782548
## 18 Epinephelus coioides 0.122874070 -0.079626206 -0.054919379
## 19 Epinephelus marginatus -0.031710533 -0.062399926 0.336729951
## 20 Euthynnus alletteratus 0.359606176 -0.123596105 -0.425146812
## 21 Gadus morhua 0.205826457 -0.231119023 -0.117400320
## 22 Genypterus blacodes -0.003802300 -0.036519512 -0.041787057
## 23 Helicolenus dactylopterus -0.372889275 0.610995446 0.512207228
## 24 Hippoglossoides platessoides -0.214739867 0.230090651 0.171685098
## 25 Hippoglossus hippoglossus -0.324693911 -0.018627963 -0.046209885
## 26 Hippoglossus stenolepis -0.060913858 0.198392074 0.126975628
## 27 Hirundichthys affinis 0.969892288 -0.681580157 -1.028682549
## 28 Hoplostethus atlanticus -0.553762785 0.740853608 0.954915957
## 29 Istiophorus platypterus 1.009444772 -0.747618418 -0.641437044
## 30 Katsuwonus pelamis 0.621543845 -0.570389152 -0.278196612
## 31 Lepidorhombus whiffiagonis -0.116472409 -0.218659779 0.008623282
## 32 Lesueurigobius friesii -0.083087339 0.155421570 0.495145392
## 33 Lethrinus lentjan 0.168784832 -0.163380349 0.086175175
## 34 Lethrinus miniatus 0.231673273 -0.054178122 0.122300241
## 35 Lethrinus nebulosus 0.099275087 -0.054486049 0.135607152
## 36 Leuresthes tenuis -0.054619350 0.074584162 0.029310216
## 37 Limanda aspera -0.214182647 0.265372484 0.380817034
## 38 Limanda ferruginea 0.232008369 -0.088670305 -0.119660388
## 39 Lithognathus mormyrus -0.100285759 -0.185985279 -0.033155000
## 40 Lopholatilus chamaeleonticeps -0.043879482 -0.165841149 0.145104671
## 41 Lutjanus analis 0.111298764 -0.214573646 0.043779287
## 42 Lutjanus carponotatus 0.198802399 -0.148932561 0.270959377
## 43 Lutjanus fulviflamma 0.028155415 0.024909304 0.323098458
## 44 Lutjanus gibbus 0.225572984 -0.117088931 0.085304142
## 45 Lutjanus purpureus -0.039061123 0.012287178 -0.137746453
## 46 Lutjanus synagris -0.008118595 -0.247260158 -0.194533328
## 47 Lutjanus vitta 0.212045030 -0.234212668 0.016125921
## 48 Maurolicus muelleri -0.004408239 -0.103867915 -0.008445044
## 49 Melanogrammus aeglefinus 0.124873982 -0.149493093 -0.033651789
## 50 Merlangius merlangus 0.134122576 -0.266939108 0.126546735
## 51 Merluccius merluccius -0.009162580 -0.144957911 -0.091193379
## 52 Micromesistius poutassou 0.002501397 -0.004141581 0.202674247
## 53 Micropogonias furnieri -0.004828232 -0.025691758 -0.432733051
## 54 Microstomus pacificus -0.223278480 0.136714569 0.578289582
## 55 Molva molva 0.211196577 -0.109859728 -0.163666752
## 56 Mullus surmuletus -0.033985423 -0.202379860 -0.005880671
## 57 Muraena helena -0.085200570 0.028312949 -0.055944934
## 58 Nemipterus japonicus 0.304769463 -0.518273151 -0.138216644
## 59 Notothenia rossii 0.072994671 0.044969527 -0.163097928
## 60 Ocyurus chrysurus 0.025818791 -0.054975119 -0.121553584
## 61 Oedalechilus labeo -0.214511226 0.250670859 0.011558889
## 62 Pagellus bogaraveo -0.171388220 0.183657103 -0.045419968
## 63 Paralabrax clathratus -0.203348260 -0.250882430 0.253675492
## 64 Pleurogrammus monopterygius 0.109008731 0.141779631 0.031499545
## 65 Pleuronectes platessa -0.182191666 -0.003357313 0.441459311
## 66 Pleuronectes quadrituberculatus -0.256610335 0.254473447 0.327448792
## 67 Pomatomus saltatrix 0.145282797 -0.256176760 -0.455572225
## 68 Pomatoschistus minutus 0.193346921 -0.270991091 -0.184565050
## 69 Pseudocyttus maculatus -0.432644454 0.736356058 0.829636070
## 70 Rachycentron canadum 0.479875160 -0.435693392 -0.370294421
## 71 Rastrelliger kanagurta 0.548204022 -0.608482335 -0.426591423
## 72 Reinhardtius hippoglossoides -0.112309040 0.238215286 0.042471535
## 73 Rhomboplites aurorubens -0.141571324 0.070715808 -0.225567340
## 74 Sardinella aurita 0.256878996 -0.109336604 -0.185453287
## 75 Sardinella gibbosa 0.486331259 -0.301602102 -0.040801296
## 76 Sardinella longiceps 0.321319798 -0.144360283 -0.397187471
## 77 Saurida tumbil -0.133819335 0.005180546 -0.469929497
## 78 Saurida undosquamis 0.352117805 -0.328875153 -0.287438508
## 79 Scomber scombrus 0.218237456 -0.238630249 0.111341937
## 80 Scomberomorus cavalla 0.199091099 -0.171265175 -0.363858765
## 81 Scomberomorus maculatus 0.325269524 -0.241252163 -0.624556196
## 82 Scophthalmus maximus 0.194822718 -0.056168604 0.127620982
## 83 Sebastes alutus -0.248880486 0.389438250 0.853016815
## 84 Sebastes caurinus -0.213855192 0.007552837 0.557033222
## 85 Sebastes crameri -0.429274239 0.574522266 0.900185509
## 86 Sebastes flavidus -0.101097569 0.265585769 0.603730530
## 87 Sebastes jordani -0.073782918 -0.089233684 0.458867350
## 88 Sebastes melanops -0.080941181 0.210475662 0.441810165
## 89 Sebastes melanostomus -0.588086024 0.549134765 0.754111016
## 90 Sebastes mentella -0.409942681 0.476249938 0.707584533
## 91 Sebastes paucispinis -0.045456509 -0.175003704 0.325549552
## 92 Sebastes polyspinis -0.182591827 0.386373401 0.705839050
## 93 Sebastes serranoides 0.081203960 0.102548199 0.310717318
## 94 Semicossyphus pulcher -0.039303533 0.013793057 0.496229911
## 95 Seriola dumerili 0.375007217 -0.267582331 -0.422732915
## 96 Sprattus sprattus 0.193951685 0.046836817 -0.029820848
## 97 Thunnus alalunga 0.239901380 -0.037294709 -0.551803703
## 98 Thunnus albacares 0.624898749 -0.663821217 -0.650559939
## 99 Thunnus maccoyii 0.331745334 -0.019298862 -0.382665143
## 100 Thunnus obesus 0.414552509 -0.283633826 -0.617115434
## 101 Thunnus orientalis 0.491071698 -0.317434186 -0.495205918
## 102 Trachinotus falcatus 0.420973532 -0.248012515 -0.037365658
## 103 Trematomus bernacchii 0.126280996 0.189634218 0.007485767
## 104 Trichiurus lepturus 0.443148474 -0.542310205 -0.340533504
## 105 Trisopterus minutus 0.075474368 -0.210073263 -0.261915471
## SpawningFreq
## 1 -0.042278790
## 2 0.675524948
## 3 0.216554864
## 4 -0.053024175
## 5 -0.058958529
## 6 -0.067040250
## 7 -0.058120576
## 8 -0.071720172
## 9 0.459697175
## 10 -0.039785515
## 11 1.566346501
## 12 -0.083653891
## 13 -0.041186509
## 14 1.215260200
## 15 0.245955627
## 16 -0.065797148
## 17 -0.055899876
## 18 -0.037256677
## 19 -0.033036514
## 20 -0.034594640
## 21 -0.025959479
## 22 -0.025617879
## 23 -0.063688760
## 24 -0.046055783
## 25 -0.003954552
## 26 -0.016572447
## 27 -0.072552278
## 28 -0.355431070
## 29 -0.002977952
## 30 -0.041082865
## 31 -0.056581428
## 32 0.182289185
## 33 -0.056860825
## 34 -0.054223751
## 35 -0.045470357
## 36 0.690768522
## 37 -0.060873930
## 38 -0.056261183
## 39 -0.064543241
## 40 -0.037101452
## 41 -0.035219643
## 42 -0.072419384
## 43 -0.071599965
## 44 -0.059025521
## 45 -0.037670573
## 46 -0.056673662
## 47 2.106829844
## 48 -0.124550432
## 49 -0.043255796
## 50 -0.058593263
## 51 -0.037746742
## 52 -0.065881739
## 53 -0.048698703
## 54 -0.059032215
## 55 -0.021529913
## 56 -0.070772387
## 57 -0.014433467
## 58 0.229686352
## 59 -0.040140775
## 60 -0.049670417
## 61 -0.076405008
## 62 -0.054090788
## 63 1.456448627
## 64 -0.061455091
## 65 -0.050644185
## 66 -0.052248579
## 67 -0.036062539
## 68 0.188960520
## 69 -0.058987649
## 70 -0.022987185
## 71 -0.072555259
## 72 -0.033684828
## 73 -0.053702454
## 74 -0.072373360
## 75 -0.086222396
## 76 -0.082127863
## 77 -0.045137501
## 78 -0.062609257
## 79 -0.063895008
## 80 -0.026472027
## 81 -0.044323836
## 82 -0.049086598
## 83 -0.059997090
## 84 -0.056261183
## 85 -0.063713407
## 86 -0.055915770
## 87 0.103903924
## 88 -0.050677776
## 89 -0.056137691
## 90 -0.059264071
## 91 0.261483055
## 92 -0.070507909
## 93 -0.059366780
## 94 1.881033301
## 95 -0.017966704
## 96 -0.093683180
## 97 -0.026849301
## 98 -0.017394342
## 99 -0.009841126
## 100 -0.012252540
## 101 -0.011028111
## 102 -0.037088994
## 103 -0.076015046
## 104 -0.025836506
## 105 -0.079043195row.names(average.resid.LHdata.phylo)<-sub(" ","_",average.resid.LHdata.phylo$Species)
average.resid.LHdata.phylo## Species K
## Anoplopoma_fimbria Anoplopoma fimbria 0.133605686
## Atractoscion_nobilis Atractoscion nobilis 0.147772321
## Aulorhynchus_flavidus Aulorhynchus flavidus 0.136152464
## Beryx_splendens Beryx splendens -0.121965077
## Caranx_crysos Caranx crysos 0.228592182
## Centropristis_striata Centropristis striata 0.250381545
## Chelidonichthys_lucerna Chelidonichthys lucerna -0.165935155
## Clupea_harengus Clupea harengus 0.060498712
## Coryphaena_hippurus Coryphaena hippurus 1.151371387
## Coryphaenoides_rupestris Coryphaenoides rupestris -0.120869106
## Cynoscion_nebulosus Cynoscion nebulosus 0.114845140
## Diplodus_annularis Diplodus annularis -0.201040215
## Engraulis_encrasicolus Engraulis encrasicolus 0.091490260
## Engraulis_mordax Engraulis mordax -0.004639927
## Eopsetta_jordani Eopsetta jordani -0.089983477
## Epinephelus_areolatus Epinephelus areolatus 0.098977047
## Epinephelus_chlorostigma Epinephelus chlorostigma 0.059012944
## Epinephelus_coioides Epinephelus coioides 0.122874070
## Epinephelus_marginatus Epinephelus marginatus -0.031710533
## Euthynnus_alletteratus Euthynnus alletteratus 0.359606176
## Gadus_morhua Gadus morhua 0.205826457
## Genypterus_blacodes Genypterus blacodes -0.003802300
## Helicolenus_dactylopterus Helicolenus dactylopterus -0.372889275
## Hippoglossoides_platessoides Hippoglossoides platessoides -0.214739867
## Hippoglossus_hippoglossus Hippoglossus hippoglossus -0.324693911
## Hippoglossus_stenolepis Hippoglossus stenolepis -0.060913858
## Hirundichthys_affinis Hirundichthys affinis 0.969892288
## Hoplostethus_atlanticus Hoplostethus atlanticus -0.553762785
## Istiophorus_platypterus Istiophorus platypterus 1.009444772
## Katsuwonus_pelamis Katsuwonus pelamis 0.621543845
## Lepidorhombus_whiffiagonis Lepidorhombus whiffiagonis -0.116472409
## Lesueurigobius_friesii Lesueurigobius friesii -0.083087339
## Lethrinus_lentjan Lethrinus lentjan 0.168784832
## Lethrinus_miniatus Lethrinus miniatus 0.231673273
## Lethrinus_nebulosus Lethrinus nebulosus 0.099275087
## Leuresthes_tenuis Leuresthes tenuis -0.054619350
## Limanda_aspera Limanda aspera -0.214182647
## Limanda_ferruginea Limanda ferruginea 0.232008369
## Lithognathus_mormyrus Lithognathus mormyrus -0.100285759
## Lopholatilus_chamaeleonticeps Lopholatilus chamaeleonticeps -0.043879482
## Lutjanus_analis Lutjanus analis 0.111298764
## Lutjanus_carponotatus Lutjanus carponotatus 0.198802399
## Lutjanus_fulviflamma Lutjanus fulviflamma 0.028155415
## Lutjanus_gibbus Lutjanus gibbus 0.225572984
## Lutjanus_purpureus Lutjanus purpureus -0.039061123
## Lutjanus_synagris Lutjanus synagris -0.008118595
## Lutjanus_vitta Lutjanus vitta 0.212045030
## Maurolicus_muelleri Maurolicus muelleri -0.004408239
## Melanogrammus_aeglefinus Melanogrammus aeglefinus 0.124873982
## Merlangius_merlangus Merlangius merlangus 0.134122576
## Merluccius_merluccius Merluccius merluccius -0.009162580
## Micromesistius_poutassou Micromesistius poutassou 0.002501397
## Micropogonias_furnieri Micropogonias furnieri -0.004828232
## Microstomus_pacificus Microstomus pacificus -0.223278480
## Molva_molva Molva molva 0.211196577
## Mullus_surmuletus Mullus surmuletus -0.033985423
## Muraena_helena Muraena helena -0.085200570
## Nemipterus_japonicus Nemipterus japonicus 0.304769463
## Notothenia_rossii Notothenia rossii 0.072994671
## Ocyurus_chrysurus Ocyurus chrysurus 0.025818791
## Oedalechilus_labeo Oedalechilus labeo -0.214511226
## Pagellus_bogaraveo Pagellus bogaraveo -0.171388220
## Paralabrax_clathratus Paralabrax clathratus -0.203348260
## Pleurogrammus_monopterygius Pleurogrammus monopterygius 0.109008731
## Pleuronectes_platessa Pleuronectes platessa -0.182191666
## Pleuronectes_quadrituberculatus Pleuronectes quadrituberculatus -0.256610335
## Pomatomus_saltatrix Pomatomus saltatrix 0.145282797
## Pomatoschistus_minutus Pomatoschistus minutus 0.193346921
## Pseudocyttus_maculatus Pseudocyttus maculatus -0.432644454
## Rachycentron_canadum Rachycentron canadum 0.479875160
## Rastrelliger_kanagurta Rastrelliger kanagurta 0.548204022
## Reinhardtius_hippoglossoides Reinhardtius hippoglossoides -0.112309040
## Rhomboplites_aurorubens Rhomboplites aurorubens -0.141571324
## Sardinella_aurita Sardinella aurita 0.256878996
## Sardinella_gibbosa Sardinella gibbosa 0.486331259
## Sardinella_longiceps Sardinella longiceps 0.321319798
## Saurida_tumbil Saurida tumbil -0.133819335
## Saurida_undosquamis Saurida undosquamis 0.352117805
## Scomber_scombrus Scomber scombrus 0.218237456
## Scomberomorus_cavalla Scomberomorus cavalla 0.199091099
## Scomberomorus_maculatus Scomberomorus maculatus 0.325269524
## Scophthalmus_maximus Scophthalmus maximus 0.194822718
## Sebastes_alutus Sebastes alutus -0.248880486
## Sebastes_caurinus Sebastes caurinus -0.213855192
## Sebastes_crameri Sebastes crameri -0.429274239
## Sebastes_flavidus Sebastes flavidus -0.101097569
## Sebastes_jordani Sebastes jordani -0.073782918
## Sebastes_melanops Sebastes melanops -0.080941181
## Sebastes_melanostomus Sebastes melanostomus -0.588086024
## Sebastes_mentella Sebastes mentella -0.409942681
## Sebastes_paucispinis Sebastes paucispinis -0.045456509
## Sebastes_polyspinis Sebastes polyspinis -0.182591827
## Sebastes_serranoides Sebastes serranoides 0.081203960
## Semicossyphus_pulcher Semicossyphus pulcher -0.039303533
## Seriola_dumerili Seriola dumerili 0.375007217
## Sprattus_sprattus Sprattus sprattus 0.193951685
## Thunnus_alalunga Thunnus alalunga 0.239901380
## Thunnus_albacares Thunnus albacares 0.624898749
## Thunnus_maccoyii Thunnus maccoyii 0.331745334
## Thunnus_obesus Thunnus obesus 0.414552509
## Thunnus_orientalis Thunnus orientalis 0.491071698
## Trachinotus_falcatus Trachinotus falcatus 0.420973532
## Trematomus_bernacchii Trematomus bernacchii 0.126280996
## Trichiurus_lepturus Trichiurus lepturus 0.443148474
## Trisopterus_minutus Trisopterus minutus 0.075474368
## tm LongevityWild SpawningFreq
## Anoplopoma_fimbria 0.036514210 0.628241328 -0.042278790
## Atractoscion_nobilis -0.228221897 -0.240573262 0.675524948
## Aulorhynchus_flavidus -0.391797789 0.050093143 0.216554864
## Beryx_splendens 0.164921460 0.129075931 -0.053024175
## Caranx_crysos -0.012072794 -0.129275401 -0.058958529
## Centropristis_striata -0.277177490 0.214774788 -0.067040250
## Chelidonichthys_lucerna -0.179376822 -0.003329178 -0.058120576
## Clupea_harengus 0.038451170 0.360566184 -0.071720172
## Coryphaena_hippurus -1.272453628 -1.002431845 0.459697175
## Coryphaenoides_rupestris 0.276869155 0.361465195 -0.039785515
## Cynoscion_nebulosus -0.286421744 -0.046178159 1.566346501
## Diplodus_annularis 0.033983409 -0.067629115 -0.083653891
## Engraulis_encrasicolus -0.251644230 -0.167626338 -0.041186509
## Engraulis_mordax -0.105730913 -0.045528672 1.215260200
## Eopsetta_jordani 0.229413680 0.332830233 0.245955627
## Epinephelus_areolatus -0.074271649 0.076851961 -0.065797148
## Epinephelus_chlorostigma 0.011153074 0.259782548 -0.055899876
## Epinephelus_coioides -0.079626206 -0.054919379 -0.037256677
## Epinephelus_marginatus -0.062399926 0.336729951 -0.033036514
## Euthynnus_alletteratus -0.123596105 -0.425146812 -0.034594640
## Gadus_morhua -0.231119023 -0.117400320 -0.025959479
## Genypterus_blacodes -0.036519512 -0.041787057 -0.025617879
## Helicolenus_dactylopterus 0.610995446 0.512207228 -0.063688760
## Hippoglossoides_platessoides 0.230090651 0.171685098 -0.046055783
## Hippoglossus_hippoglossus -0.018627963 -0.046209885 -0.003954552
## Hippoglossus_stenolepis 0.198392074 0.126975628 -0.016572447
## Hirundichthys_affinis -0.681580157 -1.028682549 -0.072552278
## Hoplostethus_atlanticus 0.740853608 0.954915957 -0.355431070
## Istiophorus_platypterus -0.747618418 -0.641437044 -0.002977952
## Katsuwonus_pelamis -0.570389152 -0.278196612 -0.041082865
## Lepidorhombus_whiffiagonis -0.218659779 0.008623282 -0.056581428
## Lesueurigobius_friesii 0.155421570 0.495145392 0.182289185
## Lethrinus_lentjan -0.163380349 0.086175175 -0.056860825
## Lethrinus_miniatus -0.054178122 0.122300241 -0.054223751
## Lethrinus_nebulosus -0.054486049 0.135607152 -0.045470357
## Leuresthes_tenuis 0.074584162 0.029310216 0.690768522
## Limanda_aspera 0.265372484 0.380817034 -0.060873930
## Limanda_ferruginea -0.088670305 -0.119660388 -0.056261183
## Lithognathus_mormyrus -0.185985279 -0.033155000 -0.064543241
## Lopholatilus_chamaeleonticeps -0.165841149 0.145104671 -0.037101452
## Lutjanus_analis -0.214573646 0.043779287 -0.035219643
## Lutjanus_carponotatus -0.148932561 0.270959377 -0.072419384
## Lutjanus_fulviflamma 0.024909304 0.323098458 -0.071599965
## Lutjanus_gibbus -0.117088931 0.085304142 -0.059025521
## Lutjanus_purpureus 0.012287178 -0.137746453 -0.037670573
## Lutjanus_synagris -0.247260158 -0.194533328 -0.056673662
## Lutjanus_vitta -0.234212668 0.016125921 2.106829844
## Maurolicus_muelleri -0.103867915 -0.008445044 -0.124550432
## Melanogrammus_aeglefinus -0.149493093 -0.033651789 -0.043255796
## Merlangius_merlangus -0.266939108 0.126546735 -0.058593263
## Merluccius_merluccius -0.144957911 -0.091193379 -0.037746742
## Micromesistius_poutassou -0.004141581 0.202674247 -0.065881739
## Micropogonias_furnieri -0.025691758 -0.432733051 -0.048698703
## Microstomus_pacificus 0.136714569 0.578289582 -0.059032215
## Molva_molva -0.109859728 -0.163666752 -0.021529913
## Mullus_surmuletus -0.202379860 -0.005880671 -0.070772387
## Muraena_helena 0.028312949 -0.055944934 -0.014433467
## Nemipterus_japonicus -0.518273151 -0.138216644 0.229686352
## Notothenia_rossii 0.044969527 -0.163097928 -0.040140775
## Ocyurus_chrysurus -0.054975119 -0.121553584 -0.049670417
## Oedalechilus_labeo 0.250670859 0.011558889 -0.076405008
## Pagellus_bogaraveo 0.183657103 -0.045419968 -0.054090788
## Paralabrax_clathratus -0.250882430 0.253675492 1.456448627
## Pleurogrammus_monopterygius 0.141779631 0.031499545 -0.061455091
## Pleuronectes_platessa -0.003357313 0.441459311 -0.050644185
## Pleuronectes_quadrituberculatus 0.254473447 0.327448792 -0.052248579
## Pomatomus_saltatrix -0.256176760 -0.455572225 -0.036062539
## Pomatoschistus_minutus -0.270991091 -0.184565050 0.188960520
## Pseudocyttus_maculatus 0.736356058 0.829636070 -0.058987649
## Rachycentron_canadum -0.435693392 -0.370294421 -0.022987185
## Rastrelliger_kanagurta -0.608482335 -0.426591423 -0.072555259
## Reinhardtius_hippoglossoides 0.238215286 0.042471535 -0.033684828
## Rhomboplites_aurorubens 0.070715808 -0.225567340 -0.053702454
## Sardinella_aurita -0.109336604 -0.185453287 -0.072373360
## Sardinella_gibbosa -0.301602102 -0.040801296 -0.086222396
## Sardinella_longiceps -0.144360283 -0.397187471 -0.082127863
## Saurida_tumbil 0.005180546 -0.469929497 -0.045137501
## Saurida_undosquamis -0.328875153 -0.287438508 -0.062609257
## Scomber_scombrus -0.238630249 0.111341937 -0.063895008
## Scomberomorus_cavalla -0.171265175 -0.363858765 -0.026472027
## Scomberomorus_maculatus -0.241252163 -0.624556196 -0.044323836
## Scophthalmus_maximus -0.056168604 0.127620982 -0.049086598
## Sebastes_alutus 0.389438250 0.853016815 -0.059997090
## Sebastes_caurinus 0.007552837 0.557033222 -0.056261183
## Sebastes_crameri 0.574522266 0.900185509 -0.063713407
## Sebastes_flavidus 0.265585769 0.603730530 -0.055915770
## Sebastes_jordani -0.089233684 0.458867350 0.103903924
## Sebastes_melanops 0.210475662 0.441810165 -0.050677776
## Sebastes_melanostomus 0.549134765 0.754111016 -0.056137691
## Sebastes_mentella 0.476249938 0.707584533 -0.059264071
## Sebastes_paucispinis -0.175003704 0.325549552 0.261483055
## Sebastes_polyspinis 0.386373401 0.705839050 -0.070507909
## Sebastes_serranoides 0.102548199 0.310717318 -0.059366780
## Semicossyphus_pulcher 0.013793057 0.496229911 1.881033301
## Seriola_dumerili -0.267582331 -0.422732915 -0.017966704
## Sprattus_sprattus 0.046836817 -0.029820848 -0.093683180
## Thunnus_alalunga -0.037294709 -0.551803703 -0.026849301
## Thunnus_albacares -0.663821217 -0.650559939 -0.017394342
## Thunnus_maccoyii -0.019298862 -0.382665143 -0.009841126
## Thunnus_obesus -0.283633826 -0.617115434 -0.012252540
## Thunnus_orientalis -0.317434186 -0.495205918 -0.011028111
## Trachinotus_falcatus -0.248012515 -0.037365658 -0.037088994
## Trematomus_bernacchii 0.189634218 0.007485767 -0.076015046
## Trichiurus_lepturus -0.542310205 -0.340533504 -0.025836506
## Trisopterus_minutus -0.210073263 -0.261915471 -0.079043195average.resid.LHdata.phylo<-average.resid.LHdata.phylo[,-1]
ppca<-phyl.pca(tree=phy,Y=average.resid.LHdata.phylo,mode="corr") #Performing the phylogenetic PCA centred and scales on species averages using the phylogeny phy
(var.expl.phylo<-round(100*diag(ppca$Eval)/sum(diag(ppca$Eval)),1))## PC1 PC2 PC3 PC4
## 56.3 25.1 13.6 5.0cumsum(var.expl.phylo) #The two first PCs explain 81.4%## PC1 PC2 PC3 PC4
## 56.3 81.4 95.0 100.0str(ppca)## List of 7
## $ Eval: num [1:4, 1:4] 2.25 0 0 0 0 ...
## ..- attr(*, "dimnames")=List of 2
## .. ..$ : chr [1:4] "PC1" "PC2" "PC3" "PC4"
## .. ..$ : chr [1:4] "PC1" "PC2" "PC3" "PC4"
## $ Evec: num [1:4, 1:4] 0.6072 -0.6047 -0.5137 0.0417 -0.0242 ...
## ..- attr(*, "dimnames")=List of 2
## .. ..$ : chr [1:4] "K" "tm" "LongevityWild" "SpawningFreq"
## .. ..$ : chr [1:4] "PC1" "PC2" "PC3" "PC4"
## $ S : num [1:105, 1:4] -0.344 0.441 0.309 -0.276 0.269 ...
## ..- attr(*, "dimnames")=List of 2
## .. ..$ : chr [1:105] "Anoplopoma_fimbria" "Atractoscion_nobilis" "Aulorhynchus_flavidus" "Beryx_splendens" ...
## .. ..$ : chr [1:4] "PC1" "PC2" "PC3" "PC4"
## $ L : num [1:4, 1:4] 0.9111 -0.9075 -0.7709 0.0626 -0.0243 ...
## ..- attr(*, "dimnames")=List of 2
## .. ..$ : chr [1:4] "K" "tm" "LongevityWild" "SpawningFreq"
## .. ..$ : chr [1:4] "PC1" "PC2" "PC3" "PC4"
## $ V : num [1:4, 1:4] 0.6309 -0.6644 -0.3244 0.0149 -0.6644 ...
## ..- attr(*, "dimnames")=List of 2
## .. ..$ : chr [1:4] "K" "tm" "LongevityWild" "SpawningFreq"
## .. ..$ : chr [1:4] "K" "tm" "LongevityWild" "SpawningFreq"
## $ a : num [1:4, 1] -0.001311 0.002492 0.000439 -0.000718
## $ mode: chr "corr"
## - attr(*, "class")= chr "phyl.pca"ppca$Evec## PC1 PC2 PC3 PC4
## K 0.60718027 -0.02423124 -0.3544356 0.71071819
## tm -0.60472681 -0.04506593 0.3741006 0.70165753
## LongevityWild -0.51371039 0.10503649 -0.8513217 0.01789913
## SpawningFreq 0.04170411 0.99315117 0.0983643 0.04728623The two first PCs explain a percentage of variance similar to the PCA without accounting for phylogeny although based here on species averages.
We then produce a biplot of the pCA by projecting single stocks on the first PC plane determined using species averages.
pcresults.phylo<-data.frame(scores(ppca,newdata=resid.LHdata.phylo[,1:4]),speciestype=resid.LHdata.phylo$speciestype) #Producing the scores for each single stock on the plane of the 2 first PCs computed using species averages.
p5<-ggplot(data=subset(pcresults.phylo,speciestype=="Other species"),aes(x=PC1,y=PC2))+
geom_point(aes(fill="Other species",shape="Other species"),color="white",size=6)+
geom_point(data=subset(pcresults.phylo,speciestype=="Tunas"),aes(fill="Tunas",shape="Tunas"),color="white",size=6)+
stat_ellipse(data=subset(pcresults.phylo,speciestype=="Tunas"),type="norm",na.rm=T,color="#CD9600",size=1)+
geom_point(data=subset(pcresults.phylo,speciestype=="Roughies"),aes(fill="Roughies",shape="Roughies"),color="white",size=6)+
stat_ellipse(data=subset(pcresults.phylo,speciestype=="Roughies"),type="norm",na.rm=T,color="#7CAE00",size=1)+
geom_segment(data=as.data.frame(ppca$Evec),aes(x=rep(0,4),y=rep(0,4),xend=5*PC1,yend=5*PC2),arrow=arrow(length=unit(0.015, "npc")),size=1)+
geom_text(data=as.data.frame(ppca$Evec),aes(x=5*PC1,y=5*PC2,label=c("Growth coefficient","Maturation age","Longevity","Spawning frequency")),nudge_x=0.75*c(1.2,-0.9,-0.6,0),nudge_y=0.2*c(0,-0.15,1,1),size=5)+
scale_shape_manual(name="Taxonomic\ngroup",values=c(
"Coelacanth circuli"=23,
"Coelacanth macrocirculi"=23,
"Coelacanth previous study"=23,
"Deep-sea sharks"=21,
"Other species"=21,
"Roughies"=21,
"Tunas"=21))+
scale_fill_manual(name="Taxonomic\ngroup",values=c(
"Coelacanth circuli"="#F8766D",
"Coelacanth macrocirculi"="#00BFC4",
"Coelacanth previous study"="#619CFF",
"Deep-sea sharks"="#C77CFF",
"Other species"="darkgrey",
"Roughies"="#7CAE00",
"Tunas"="#CD9600"))+
xlab(paste("PC1 (",var.expl.phylo[1],"%)",sep=""))+
ylab(paste("PC2 (",var.expl.phylo[2],"%)",sep=""))+
scale_x_continuous(breaks=c(-4,-2,0,2,4))+
scale_y_continuous(breaks=c(-2.5,0,2.5,5,7.5))+
theme(legend.position="left",
legend.key=element_blank(),
legend.title=element_text(size=14),
legend.text=element_text(size=12),
axis.title.x = element_text(size=16, face="bold"),
axis.title.y = element_text(size=16, face="bold"),
axis.text=element_text(size=12))+
geom_hline(yintercept=0,linetype=2)+
geom_vline(xintercept=0,linetype=2)+
guides(shape = guide_legend(override.aes = list(size = 5)))As for the analysis without phylogeny, we produce a boxplot of the cores along the first PC which is an index of the slow-fast continuum.
p6<-ggplot(data=pcresults.phylo, aes(x=speciestype, y=PC1)) +
geom_boxplot(aes(fill=speciestype),show.legend=F,outlier.size=3)+
stat_summary(fun=mean, geom="point", shape=21, fill="black",size=4) +
scale_fill_manual(values=c(
"Deep-sea sharks"="#C77CFF",
"Other species"="darkgrey",
"Roughies"="#7CAE00",
"Tunas"="#CD9600"))+
labs(x="Taxonomic group", y="Slow-fast life-history axis, PC1") +
coord_flip()+
theme(
axis.title.x = element_text(size=16, face="bold"),
axis.title.y = element_text(size=16, face="bold"),
axis.text=element_text(size=12),
legend.key = element_blank()
)Then we combine the plots.
ggarrange(p5,p6,labels=c("A","B"),ncol=1,nrow=2,heights=c(2,1),common.legend=T,legend="left") 
As we can see, results are very similar to the PCA without accounting for the phylogeny.