knitr::opts_chunk$set(collapse = TRUE, comment = ">", dev = 'pdf')
RClone data format: several populations
RClone is a R package version of GenClone program (Arnaud-Haond & Belkhir 2007): to analyse data (SSR, SNP, ...), test for clonality and describe spatial clonal organisation.
Major improvements are multi-populations handling and definition of MLLs (Multilocus Lineages, i.e. slightly distinct Multi Locus Genotypes) through simulations.
RClone allows:
Some of these analysis can be applied to dataset with no repeated MLG, regardless of the reproductive system (sexual, partially asexual or strictly asexual).
RClone functions works on diploid/haploid, one or several populations dataset.
If you have only one population in your dataset, go to other vignette RClone_quickmanual.
If you have haploid data, you can skip to D. Description of data set.
To use RClone functions, your data table must look like:
library(RClone) data(posidonia)
knitr::kable(posidonia[1:10,1:8], align = "c")
There is only one allele per column and, per locus, alleles are sorted by increasing order.
This is mandatory for all RClone functions.
As formatting can be source of error, we included functions to help formatting your diploid data:
1, The classic infile you could have, one locus per column
data(zostera) head(zostera)
data(zostera) knitr::kable(head(zostera), align = "c")
Zostera data is composed of:
a first column with population indication;
a second and third columns with x/y coordinates;
* a genotypic dataset.
popvec <- zostera[,1] #futur vecpop coord_zostera <- zostera[,2:3] #futur coordinates zostera <- zostera[,4:ncol(zostera)] #dataset zostera <- convert_GC(zostera, 3) #We used "3" because this is the length of each allele.
head(zostera)
knitr::kable(zostera[1:6,1:7], align = "c")
2, The simple case: you already have a one-allele per column table
Just remove the pop/coords informations as above and sort your alleles:
sort_all(zostera)
3, You already work with Adegenet
Similar to case number 1, except you have to export your genind
data into table first:
#library(adegenet) #with data1, a genind object from Adegenet: test <- genind2df(data1) data2 <- convert_GC(test, 3, "/") #only if yours alleles are of length "3"
List unique alleles per locus:
Basic commands:
list_all_tab(zostera, vecpop = popvec)
or, for haploid data:
list_all_tab(haplodata, haploid = TRUE, vecpop = haplovec)
Results:
list_all_tab(zostera, vecpop = popvec)
#SaintMalo
knitr::kable(list_all_tab(zostera, vecpop = popvec)[[1]], align = "c")
#Arcouest
knitr::kable(list_all_tab(zostera, vecpop = popvec)[[2]], align = "c")
List MLG:
Basic commands:
MLG_tab(zostera, vecpop = popvec)
or, for haploid data:
MLG_tab(haplodata, vecpop = haplovec)
Results:
MLG_tab(zostera, vecpop = popvec)[[1]] #SaintMalo
knitr::kable(MLG_tab(zostera, vecpop = popvec)[[1]][1:5,], align = "c")
Allelic frequencies:
Basic commands:
freq_RR(zostera, vecpop = popvec)
or, for haploid data:
freq_RR(haplodata, haploid = TRUE, vecpop = haplovec)
Options:
freq_RR(zostera, vecpop = popvec) #on ramets freq_RR(zostera, vecpop = popvec, genet = TRUE) #on genets freq_RR(zostera, vecpop = popvec, RR = TRUE) #Round-Robin methods
Results:
freq_RR(zostera, vecpop = popvec)[[1]] #SaintMalo
res <- cbind(freq_RR(zostera, vecpop = popvec)[[1]], freq_RR(zostera, vecpop = popvec, genet = TRUE)[[1]][,3], freq_RR(zostera, vecpop = popvec, RR = TRUE)[[1]][,3])[1:7,] colnames(res)[3:5] <- c("freq_ramet", "freq_genet", "freq_RR") knitr::kable(res, align = "c")
On loci
Basic commands:
sample_loci(zostera, vecpop = popvec, nbrepeat = 1000)
or, for haploid data:
sample_loci(haplodata, haploid = TRUE, vecpop = haplovec, nbrepeat = 1000)
Options:
sample_loci(zostera, vecpop = popvec, nbrepeat = 1000, He = TRUE) #with He results sample_loci(zostera, vecpop = popvec, nbrepeat = 1000, graph = TRUE) #graph displayed sample_loci(zostera, vecpop = popvec, nbrepeat = 1000, bar = TRUE) #progression bar, could be time consuming sample_loci(zostera, vecpop = popvec, nbrepeat = 1000, export = TRUE) #graph export in .eps
Results:
res <- sample_loci(zostera, vecpop = popvec, nbrepeat = 1000, He = TRUE) names(res)
data(resvigncont2) names(resvigncont2$res2_SU1)
names(res$SaintMalo)
names(resvigncont2$res2_SU1$SaintMalo)
#Results: MLG res$Arcouest$res_MLG
knitr::kable(resvigncont2$res2_SU1$Arcouest$res_MLG, align = "c")
#Results: alleles res$Arcouest$res_alleles
knitr::kable(res$Arcouest$res_alleles, align = "c")
#Results: raw data #res$Arcouest$raw_He #res$Arcouest$raw_MLG #res$Arcouest$raw_all
boxplot(res$SaintMalo$raw_MLG, main = "Genotype accumulation curve", xlab = "Number of loci sampled", ylab = "Number of multilocus genotypes")
boxplot(resvigncont2$res2_SU1$SaintMalo$raw_MLG, main = "Genotype accumulation curve", xlab = "Number of loci sampled", ylab = "Number of multilocus genotypes")
boxplot(res$Arcouest$raw_MLG, main = "Genotype accumulation curve", xlab = "Number of loci sampled", ylab = "Number of multilocus genotypes")
boxplot(resvigncont2$res2_SU1$Arcouest$raw_MLG, main = "Genotype accumulation curve", xlab = "Number of loci sampled", ylab = "Number of multilocus genotypes")
Same on units
Basic commands:
sample_units(zostera, vecpop = popvec, nbrepeat = 1000)
or, for haploid data:
sample_units(haplodata, haploid = TRUE, vecpop = haplovec, nbrepeat = 1000)
This sub-sampling analysis deliver basic estimates of richness and diversity for an increasing number of sampling units.
They can be used to standardise estimates of populations with different sampling effort.
pgen, psex and p-values
Basic commands:
pgen(zostera, vecpop = popvec) psex(zostera, vecpop = popvec)
or, for haploid data:
pgen(haplodata, haploid = TRUE, vecpop = haplovec) psex(haplodata, haploid = TRUE, vecpop = haplovec)
Options: (idem on psex and pgen)
#allelic frequencies computation: psex(zostera, vecpop = popvec) #psex on ramets psex(zostera, vecpop = popvec, genet = TRUE) #psex on genets psex(zostera, vecpop = popvec, RR = TRUE) #psex with Round-Robin method #psex computation psex(zostera, vecpop = popvec) #psex with one psex per replica psex(zostera, vecpop = popvec, MLGsim = TRUE) #psex MLGsim method #pvalues: psex(zostera, vecpop = popvec, nbrepeat = 100) #with p-values psex(zostera, vecpop = popvec, nbrepeat = 1000, bar = TRUE) #with p-values and a progression bar
Results:
res <- psex(zostera, vecpop = popvec, RR = TRUE, nbrepeat = 1000) res$Arcouest[[1]] #if nbrepeat != 0, res contains a table of psex values and a vector of sim-psex values
knitr::kable(resvigncont2$res2_PS2, align = "c")
res$Arcouest[[2]] #a part of sim-psex values
resvigncont2$res2_PS1$Arcouest[[2]][1:10]
Fis, pgen Fis, psex Fis and p-values
Not for haploid data !
Fis
Basic commands:
Fis(zostera, vecpop = popvec)
Options:
Fis(zostera, vecpop = popvec) #Fis on ramets Fis(zostera, vecpop = popvec, genet = TRUE) #Fis on genets Fis(zostera, vecpop = popvec, RR = TRUE) #Fis with Round-Robin methods #RR = TRUE contains two results : a table with allelic frequencies #and a table with Fis results
Results:
Fis(zostera, vecpop = popvec, RR = TRUE)$Arcouest[[2]]
knitr::kable(Fis(zostera, vecpop = popvec, RR = TRUE)$Arcouest[[2]], align = "c")
pgen Fis, psex Fis and p-values
Basic commands: (idem for pgen_Fis and psex_Fis)
pgen_Fis(zostera, vecpop = popvec)
Options:
#allelic frequencies: psex_Fis(zostera, vecpop = popvec) #psex Fis on ramets psex_Fis(zostera, vecpop = popvec, genet = TRUE) #psex Fis on genets psex_Fis(zostera, vecpop = popvec, RR = TRUE) #psex Fis with Round-Robin method #psex computation psex_Fis(zostera, vecpop = popvec) #psex Fis, one for each replica psex_Fis(zostera, vecpop = popvec, MLGsim = TRUE) #psex Fis with MLGsim method #pvalues psex_Fis(zostera, vecpop = popvec, nbrepeat = 100) #with p-values psex_Fis(zostera, vecpop = popvec, nbrepeat = 1000, bar = TRUE) #with p-values and a progression bar
Results:
res <- psex_Fis(zostera, vecpop = popvec, RR = TRUE, nbrepeat = 1000) res$Arcouest[[1]] #if nbrepeat != 0, res contains a table of psex values and a vector of sim-psex Fis values
knitr::kable(resvigncont2$res2_PS4, align = "c")
res$Arcouest[[2]] #a part of sim psex Fis values
resvigncont2$res2_PS3$Arcouest[[2]][1:10]
Genetic distance matrix computation and threshold definition
On a theoretical diploid population with c = 0.9999 (c, clonality rate):
data(popsim) vecsim <- c(rep(1,50), rep(2,50))
#genetic distances computation, distance on allele differences: respop <- genet_dist(popsim, vecpop = vecsim) ressim <- genet_dist_sim(popsim, vecpop = vecsim , nbrepeat = 1000) #theoretical distribution: sexual reproduction ressimWS <- genet_dist_sim(popsim, vecpop = vecsim , genet = TRUE, nbrepeat = 1000) #idem, without selfing
respop <- resvigncont2$respop ressim <- resvigncont2$ressim ressimWS <- resvigncont2$ressimWS
#graph prep.: #first pop: p1 <- hist(respop[[1]]$distance_matrix, freq = FALSE, col = rgb(0,0.4,1,1), breaks = seq(0, max(respop[[1]]$distance_matrix)+1, 1), main = "pop_1_sim", xlab = "") p2 <- hist(ressim[[1]]$distance_matrix, freq = FALSE, col = rgb(0.7,0.9,1,0.5), breaks = seq(0, max(ressim[[1]]$distance_matrix)+1, 1), main = "pop_1_SR", xlab = "") p3 <- hist(ressimWS[[1]]$distance_matrix, freq = FALSE, col = rgb(0.9,0.5,1,0.3), breaks = seq(0, max(ressimWS[[1]]$distance_matrix)+1, 1), main = "pop_1_SRWS", xlab = "") limx <- max(max(respop[[1]]$distance_matrix), max(ressim[[1]]$distance_matrix), max(ressimWS[[1]]$distance_matrix)) #graph superposition: plot(p1, col = rgb(0,0.4,1,1), freq = FALSE, xlim = c(0,limx), main = paste("pop", unique(vecsim)[[1]], sep = "_"), xlab = "Genetic distances") plot(p2, col = rgb(0.7,0.9,1,0.5), freq = FALSE, add = TRUE) plot(p3, col = rgb(0.9,0.5,1,0.3), freq = FALSE, add = TRUE) #adding a legend: leg.txt <- c("original data","simulated data", "without selfing") col <- c(rgb(0,0.4,1,1), rgb(0.7,0.9,1,0.5), rgb(0.9,0.5,1,0.3)) legend("top", fill = col, leg.txt, plot = TRUE, bty = "o", box.lwd = 1.5, bg = "white") #second pop: p <- 2 #useful if several populations: just change *p* and run lines p1 <- hist(respop[[p]]$distance_matrix, freq = FALSE, col = rgb(0,0.4,1,1), breaks = seq(0, max(respop[[p]]$distance_matrix)+1, 1), main = paste("pop", p, sep = "_"), xlab = "") p2 <- hist(ressim[[p]]$distance_matrix, freq = FALSE, col = rgb(0.7,0.9,1,0.5), breaks = seq(0, max(ressim[[p]]$distance_matrix)+1, 1), main = paste("pop", p, sep = "_"), xlab = "") p3 <- hist(ressimWS[[p]]$distance_matrix, freq = FALSE, col = rgb(0.9,0.5,1,0.3), breaks = seq(0, max(ressimWS[[p]]$distance_matrix)+1, 1), main = paste("pop", p, sep = "_"), xlab = "") limx <- max(max(respop[[p]]$distance_matrix), max(ressim[[p]]$distance_matrix), max(ressimWS[[p]]$distance_matrix)) #graph superposition: plot(p1, col = rgb(0,0.4,1,1), freq = FALSE, xlim = c(0,limx), main = paste("pop", unique(vecsim)[[p]], sep = "_"), xlab = "Genetic distances") plot(p2, col = rgb(0.7,0.9,1,0.5), freq = FALSE, add = TRUE) plot(p3, col = rgb(0.9,0.5,1,0.3), freq = FALSE, add = TRUE) #adding a legend: leg.txt <- c("original data","simulated data", "without selfing") col <- c(rgb(0,0.4,1,1), rgb(0.7,0.9,1,0.5), rgb(0.9,0.5,1,0.3)) legend("top", fill = col, leg.txt, plot = TRUE, bty = "o", box.lwd = 1.5, bg = "white")
#determining alpha2 table(respop[[1]]$distance_matrix) #alpha2 = 3
#creating MLL list: MLLlist <- MLL_generator(popsim, vecpop = vecsim, alpha2 = c(3,0)) ##This will create a list of MLL (alpha2 = 3) and MLG (alpha2 = 0) ! #or res <- genet_dist(popsim, vecpop = vecsim, alpha2 = c(3,0)) MLLlist <- MLL_generator2(list(res[[1]]$potential_clones, res[[2]]$potential_clones), MLG_list(popsim, vecpop = vecsim), vecpop = vecsim)
For haploid data, theoretical example:
respop <- genet_dist(haplodata, haploid = TRUE, vecpop = vechaplo) ressim <- genet_dist_sim(haplodata, haploid = TRUE, vecpop = vechaplo, nbrepeat = 1000) MLLlist <- MLL_generator(haplodata, haploid = TRUE, vecpop = vechaplo, alpha2 = c(3,0)) #or res <- genet_dist(haplodata, haploid = TRUE, vecpop = vechaplo, alpha2 = c(3,0)) MLLlist <- MLL_generator2(list(res[[1]]$potential_clones, res[[2]]$potential_clones), haploid = TRUE, MLG_list(haplodata, vecpop = vechaplo), vecpop = vechaplo)
Basic commands:
clonal_index(zostera, vecpop = popvec)
or, with MLL:
clonal_index(popsim, vecpop = vecsim, listMLL = MLLlist)
or, for haploid data:
clonal_index(haplodata, vecpop = vechaplo)
Results:
clonal_index(zostera, vecpop = popvec)
knitr::kable(resvigncont2$res2_ci, align = "c")
Basic commands:
Pareto_index(zostera, vecpop = popvec)
or, with MLL:
Pareto_index(popsim, vecpop = vecsim, listMLL = MLLlist)
or, for haploid data:
Pareto_index(haplodata, vecpop = vechaplo)
Options:
Pareto_index(zostera, vecpop = popvec, graph = TRUE) #classic graphic Pareto_index(zostera, vecpop = popvec, legends = 2, export = TRUE) #export option Pareto_index(zostera, vecpop = popvec, full = TRUE) #all results
Results:
res <- Pareto_index(zostera, vecpop = popvec, full = TRUE, graph = TRUE, legends = 2)
require(RClone) resz <- split(zostera, popvec) res1 <- Pareto_index(resz[[1]], full = TRUE, graph = TRUE, legends = 2) res2 <- Pareto_index(resz[[2]], full = TRUE, graph = TRUE, legends = 2)
names(res$SaintMalo)
names(res1)
res$SaintMalo$Pareto
res1$Pareto
res$SaintMalo$c_Pareto
res1$c_Pareto
#res$SaintMalo$regression_results #res$SaintMalo$coords_Pareto #points coordinates
Basic commands:
autocorrelation(zostera, coords = coord_zostera, vecpop = popvec, Loiselle = TRUE)
or, with MLL:
autocorrelation(popsim, coords = coord_sim, Loiselle = TRUE, listMLL = MLLlist)
or, for haploid data:
autocorrelation(haplodata, haploid = TRUE, coords = coord_haplo, Loiselle = TRUE)
Lot's of options:
#kinship distances: autocorrelation(zostera, coords = coord_zostera, vecpop = popvec, Loiselle = TRUE) autocorrelation(zostera, coords = coord_zostera, vecpop = popvec, Ritland = TRUE) #ramets/genets methods: autocorrelation(zostera, coords = coord_zostera, vecpop = popvec, Loiselle = TRUE) #ramets autocorrelation(zostera, coords = coord_zostera, vecpop = popvec,, Loiselle = TRUE, genet = TRUE, central_coords = TRUE) #genets, central coordinates of each MLG autocorrelation(zostera, coords = coord_zostera, vecpop = popvec, Loiselle = TRUE, genet = TRUE, random_unit = TRUE) #genets, one random unit per MLG autocorrelation(zostera, coords = coord_zostera, vecpop = popvec, Loiselle = TRUE, genet = TRUE, weighted = TRUE) #genets, with weighted matrix on kinships #distance classes construction: autocorrelation(zostera, coords = coord_zostera, vecpop = popvec, Loiselle = TRUE) #10 equidistant classes distvec <- c(0,10,15,20,30,50,70,76.0411074) #with 0, min distance and 76.0411074, max distance autocorrelation(zostera, coords = coord_zostera, vecpop = popvec, Loiselle = TRUE, vecdist = distvec) #custom distance vector autocorrelation(zostera, coords = coord_zostera, vecpop = popvec, Loiselle = TRUE, class1 = TRUE, d = 7) #7 equidistant classes autocorrelation(zostera, coords = coord_zostera, vecpop = popvec, Loiselle = TRUE, class2 = TRUE, d = 7) #7 distance classes with the same number of units in each #graph options: autocorrelation(zostera, coords = coord_zostera, vecpop = popvec, Ritland = TRUE, graph = TRUE) #displays graph autocorrelation(zostera, coords = coord_zostera, vecpop = popvec, Ritland = TRUE, export = TRUE) #export graph #pvalues computation autocorrelation(zostera, coords = coord_zostera, vecpop = popvec, Ritland = TRUE, nbrepeat = 1000)
Results:
res <- autocorrelation(zostera, coords = coord_zostera, vecpop = popvec, Ritland = TRUE, nbrepeat = 1000, graph = TRUE)
plot(resvigncont2$res2auto1$Main_results[,3], resvigncont2$res2auto1$Main_results[,6], main = "Spatial aucorrelation analysis", ylim = c(-0.2,0.2), type = "l", xlab = "Spatial distance", ylab = "Coancestry (Fij)") points(resvigncont2$res2auto1$Main_results[,3], resvigncont2$res2auto1$Main_results[,6], pch = 20) abline(h = 0, lty = 3) plot(resvigncont2$res2auto2$Main_results[,3], resvigncont2$res2auto2$Main_results[,6], main = "Spatial aucorrelation analysis", ylim = c(-0.2,0.2), type = "l", xlab = "Spatial distance", ylab = "Coancestry (Fij)") points(resvigncont2$res2auto2$Main_results[,3], resvigncont2$res2auto2$Main_results[,6], pch = 20) abline(h = 0, lty = 3)
names(res$Arcouest)
names(resvigncont2$res2auto2)
res$Arcouest$Main_results #enables graph reproduction
knitr::kable(resvigncont2$res2auto2$Main_results, align = "c")
apply(res$Arcouest$Main_results, 2, mean)[6] #mean Fij
apply(resvigncont2$res2auto2$Main_results, 2, mean)[6]
res$Arcouest$Slope_and_Sp_index #gives b and Sp indices
knitr::kable(resvigncont2$res2auto2$Slope_and_Sp_index, align = "c")
#raw data: #res$Arcouest$Slope_resample #res$Arcouest$Kinship_resample #res$Arcouest$Matrix_kinship_results #res$Arcouest$Class_kinship_results #res$Arcouest$Class_distance_results
Basic commands:
clonal_sub(zostera, coords = coord_zostera, vecpop = popvec)
or, with MLL:
clonal_sub(popsim, coords = coord_sim, listMLL = MLLlist)
or, for haploid data:
clonal_sub(haplodata, haploid = TRUE, coords = coord_haplo)
Options: same distance classes definition as autocorrelation:
clonal_sub(posidonia, coords = coord_posidonia) #basic, with 10 equidistant classes distvec <- c(0,10,15,20,30,50,70,76.0411074) #with 0, min distance and 76.0411074, max distance clonal_sub(zostera, coords = coord_zostera, vecpop = popvec, vecdist = distvec) #custom distance classes clonal_sub(zostera, coords = coord_zostera, vecpop = popvec, class1 = TRUE, d = 7) #7 equidistant classes clonal_sub(zostera, coords = coord_zostera, vecpop = popvec, class1 = TRUE, d = 7) #7 distance classes with the same number of units in each
Results:
res <- clonal_sub(zostera, coords = coord_zostera, vecpop = popvec) res$Arcouest[[1]] #Global clonal subrange
res$Arcouest$clonal_sub_tab #details per class
knitr::kable(res$Arcouest$clonal_sub_tab, align ="c")
Basic commands:
agg_index(zostera, coords = coord_zostera, vecpop = popvec)
or, with MLL:
agg_index(popsim, coords = coord_sim, listMLL = MLLlist)
or, for haploid data:
agg_index(haplodata, coords = coord_haplo)
Options:
agg_index(zostera, coords = coord_zostera, vecpop = popvec, nbrepeat = 100) #pvalue computation agg_index(zostera, coords = coord_zostera, vecpop = popvec, nbrepeat = 1000, bar = TRUE) #could be time consuming
Results:
res <- agg_index(zostera, coords = coord_zostera, vecpop = popvec, nbrepeat = 1000)
res$SaintMalo$results #Aggregation index
knitr::kable(resvigncont2$res2_agg$SaintMalo$results, align = "c")
#res$SaintMalo$simulation #vector of sim aggregation index
Basic commands:
#for zostera, centers of quadra is at 15,10 edge_effect(zostera, coords = coord_zostera, vecpop = popvec, center = rep(c(15,10),2))
or, with MLL:
edge_effect(popsim, coords = coord_sim, center = rep(c(15,10),2), listMLL = MLLlist)
or, for haploid data:
edge_effect(haplodata, coords = coord_haplo, center = rep(c(15,10),2))
Options:
edge_effect(zostera, coords = coord_zostera, vecpop = popvec, center = rep(c(15,10),2), nbrepeat = 100) #pvalue computation edge_effect(zostera, coords = coord_zostera, vecpop = popvec, center = rep(c(15,10),2), nbrepeat = 1000, bar = TRUE) #could be time consuming
Results:
res <- edge_effect(zostera, coords = coord_zostera, vecpop = popvec, center = rep(c(15,10),2), nbrepeat = 100) #better put 1000 nbrepeat at least
res$SaintMalo$results #Aggregation index
knitr::kable(resvigncont2$res2_ee$SaintMalo$results, align = "c")
#res$SaintMalo$simulation #vector of sim aggregation index
Summary function of main results:
Basic commands:
GenClone(zostera, coords = coord_zostera, vecpop = popvec)
or, with MLL:
GenClone(popsim, coords = coord_sim, listMLL = MLLlist)
or, for haploid data:
GenClone(haplodata, haploid = TRUE, coords = coord_haplo)
Options:
GenClone(zostera, coords = coord_zostera, vecpop = popvec, nbrepeat = 100) #pvalues GenClone(zostera, coords = coord_zostera, vecpop = popvec, nbrepeat = 1000, bar = TRUE) #could be time consuming
Results:
GenClone(zostera, coords = coord_zostera, vecpop = popvec)
knitr::kable(resvigncont2$res2_gen[,1:9], longtable = TRUE, align = "c") knitr::kable(resvigncont2$res2_gen[,10:16], longtable = TRUE, align = "c") knitr::kable(resvigncont2$res2_gen[,17:24], longtable = TRUE, align = "c")
When a locus is homozygous, it is ignored for Fis and Fis_WR values computation.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.