Hello everyone, I am having a lot of trouble doing the RDA chart. I'm new in programming and R. I don't really know if I'm doing things right. I am trying to associate certain parameters with the presence of microorganisms, but I want this to be associated according to the microorganisms found in July, August and September. I enclose my script and the image I obtained so far. I don't know how to make Julio's microorganisms paint orange and so on. 
# Sample data on a copy/paste friendly format
# clear my environment
rm(list=ls())
#RDA
rm(list = ls())
metadata <- data.frame(stringsAsFactors=FALSE,
SampleID = c("1A", "1B", "1C", "1D", "1E", "1F", "1G", "1H", "2A",
"2B", "2C", "2D", "2E", "2F", "2G", "2H", "3A", "3B",
"3C", "3D", "3E", "3F", "3G", "3H", "4A", "4B", "4C",
"4D", "4E", "4F", "4G", "4H", "5A", "5B", "5C", "5D",
"5E", "5F", "5G", "5H", "6A", "6B", "6C", "6D", "6E", "6F",
"6G", "6H", "7A", "7B", "7C", "7D", "7E", "7F", "7G",
"7H", "8A", "8B", "8C", "8D", "8E", "8F", "8G", "8H",
"9A", "9B", "9C", "9D", "9E", "9F", "9G", "9H", "10A",
"10B", "10C", "10D", "10E", "10F", "10G", "10H", "11A",
"11B", "11C", "11D", "11E", "11F", "11G", "11H", "12A",
"12B", "12C", "12D", "12E", "12F", "12G", "12H"),
SamplingPoint = c("CAJ-1", "CAJ-2", "CAJ-4", "CAJ-1", "CAJ-3", "CAJ-5",
"CAJ-2", "CAJ-4", "CAJ-1", "CAJ-2", "CAJ-4", "CAJ-1",
"CAJ-3", "CAJ-5", "CAJ-2", "CAJ-4", "CAJ-1", "CAJ-2",
"CAJ-4", "CAJ-1", "CAJ-3", "CAJ-5", "CAJ-2", "CAJ-4",
"CAJ-1", "CAJ-3", "CAJ-4", "CAJ-1", "CAJ-3", "CAJ-5",
"CAJ-2", "CAJ-4", "CAJ-1", "CAJ-3", "CAJ-4", "CAJ-1", "CAJ-3",
"CAJ-5", "CAJ-2", "CAJ-4", "CAJ-1", "CAJ-3", "CAJ-5",
"CAJ-1", "CAJ-3", "CAJ-5", "CAJ-2", "CAJ-4", "CAJ-1",
"CAJ-3", "CAJ-5", "CAJ-2", "CAJ-4", "CAJ-1", "CAJ-3",
"CAJ-5", "CAJ-1", "CAJ-3", "CAJ-5", "CAJ-2", "CAJ-4",
"CAJ-1", "CAJ-3", "CAJ-5", "CAJ-2", "CAJ-3", "CAJ-5", "CAJ-2",
"CAJ-4", "CAJ-1", "CAJ-3", "CAJ-5", "CAJ-2", "CAJ-3",
"CAJ-5", "CAJ-2", "CAJ-4", "CAJ-1", "CAJ-3", "CAJ-5",
"CAJ-2", "CAJ-4", "CAJ-5", "CAJ-2", "CAJ-4", "CAJ-1",
"CAJ-3", "CAJ-5", "CAJ-2", "CAJ-4", "CAJ-1", "CAJ-2", "CAJ-4",
"CAJ-1", "CAJ-3", "CAJ-5"),
Depth = c("80 cm", "80 cm", "Interstitial", "80 cm", "80 cm",
"80 cm", "80 cm", "80 cm", "80 cm", "Interstitial",
"Interstitial", "80 cm", "80 cm", "80 cm", "80 cm",
"80 cm", "Interstitial", "Interstitial", "80 cm",
"Interstitial", "Interstitial", "Interstitial", "Interstitial",
"Interstitial", "Interstitial", "80 cm", "Interstitial",
"Interstitial", "Interstitial", "Interstitial",
"Interstitial", "Interstitial", "80 cm", "80 cm", "Interstitial",
"80 cm", "80 cm", "80 cm", "80 cm", "80 cm", "80 cm",
"Interstitial", "80 cm", "80 cm", "Interstitial",
"Interstitial", "Interstitial", "Interstitial", "Interstitial",
"Interstitial", "80 cm", "80 cm", "80 cm", "80 cm",
"80 cm", "80 cm", "Interstitial", "80 cm", "Interstitial",
"80 cm", "80 cm", "80 cm", "80 cm", "80 cm", "80 cm",
"Interstitial", "Interstitial", "Interstitial",
"Interstitial", "Interstitial", "Interstitial", "Interstitial",
"80 cm", "Interstitial", "80 cm", "Interstitial",
"Interstitial", "Interstitial", "Interstitial", "Interstitial",
"Interstitial", "80 cm", "Interstitial", "80 cm", "80 cm",
"80 cm", "80 cm", "80 cm", "Interstitial", "80 cm",
"Interstitial", "Interstitial", "Interstitial",
"Interstitial", "Interstitial", "Interstitial"),
Month = c("July", "July", "July", "August", "August", "August",
"September", "September", "July", "July", "July",
"August", "August", "August", "September", "September",
"July", "July", "July", "August", "August", "August",
"September", "September", "July", "July", "July", "August",
"August", "August", "September", "September", "July",
"July", "July", "August", "August", "August", "September",
"September", "July", "July", "July", "August", "August",
"August", "September", "September", "July", "July",
"July", "August", "August", "September", "September",
"September", "July", "July", "July", "August", "August",
"September", "September", "September", "July", "July",
"July", "August", "August", "September", "September",
"September", "July", "July", "July", "August", "August",
"September", "September", "September", "July", "July", "July",
"August", "August", "September", "September",
"September", "July", "July", "August", "August", "August",
"September", "September", "September"),
Filter = c("20-25 um", "0.45 um", "20-25 um", "20-25 um",
"20-25 um", "20-25 um", "20-25 um", "20-25 um",
"0.45 um", "20-25 um", "0.45 um", "0.45 um", "0.45 um",
"0.45 um", "0.45 um", "0.45 um", "20-25 um", "0.45 um",
"0.45 um", "20-25 um", "20-25 um", "20-25 um", "20-25 um",
"20-25 um", "0.45 um", "20-25 um", "20-25 um", "0.45 um",
"0.45 um", "0.45 um", "0.45 um", "0.45 um", "20-25 um",
"0.45 um", "0.45 um", "0.45 um", "0.45 um", "0.45 um",
"0.45 um", "0.45 um", "0.45 um", "20-25 um", "20-25 um",
"20-25 um", "0.45 um", "0.45 um", "0.45 um", "0.45 um",
"20-25 um", "0.45 um", "0.45 um", "20-25 um", "20-25 um",
"20-25 um", "20-25 um", "20-25 um", "0.45 um",
"0.45 um", "20-25 um", "0.45 um", "0.45 um", "0.45 um",
"0.45 um", "0.45 um", "20-25 um", "20-25 um", "0.45 um",
"20-25 um", "20-25 um", "20-25 um", "20-25 um", "20-25 um",
"0.45 um", "0.45 um", "0.45 um", "0.45 um", "0.45 um",
"0.45 um", "0.45 um", "0.45 um", "20-25 um", "20-25 um",
"20-25 um", "0.45 um", "0.45 um", "0.45 um", "0.45 um",
"0.45 um", "0.45 um", "0.45 um", "0.45 um", "0.45 um",
"0.45 um", "0.45 um", "0.45 um", "0.45 um"),
SampleorReplica = c("Sample", "Sample", "Sample", "Sample", "Sample",
"Sample", "Sample", "Sample", "Sample", "Sample",
"Sample", "Sample", "Sample", "Sample", "Sample", "Sample",
"Sample", "Sample", "Sample", "Sample", "Sample",
"Sample", "Sample", "Sample", "Sample", "Sample", "Replica",
"Sample", "Sample", "Sample", "Sample", "Sample",
"Replica", "Sample", "Replica", "Replica", "Replica",
"Replica", "Replica", "Replica", "Replica", "Sample", "Sample",
"Replica", "Replica", "Replica", "Replica", "Replica",
"Replica", "Sample", "Sample", "Sample", "Sample",
"Sample", "Sample", "Sample", "Replica", "Replica", "Sample",
"Sample", "Sample", "Sample", "Sample", "Sample",
"Sample", "Replica", "Sample", "Sample", "Sample", "Sample",
"Sample", "Sample", "Replica", "Replica", "Replica",
"Sample", "Sample", "Sample", "Sample", "Sample", "Replica",
"Sample", "Replica", "Replica", "Replica", "Replica",
"Replica", "Replica", "Replica", "Replica", "Replica",
"Replica", "Replica", "Replica", "Replica", "Replica"),
DO = c(4.24, 2.58, 2.85, 3.81, 2.64, 4.61, 5.74, 3.69, 4.24,
2.66, 2.85, 3.81, 2.64, 4.61, 5.74, 3.69, 2.51, 2.66,
3.14, 2.92, 1.83, 1.82, 3.02, 3.3, 2.51, 3.04, 2.85,
2.92, 1.83, 1.82, 3.02, 3.3, 4.24, 3.04, 2.85, 3.81, 2.64,
4.61, 5.74, 3.69, 4.24, 2.09, 2.98, 3.81, 1.83, 1.82,
3.02, 3.3, 2.51, 2.09, 2.98, 8.17, 2.63, 7.69, 3.07, 3.7,
2.51, 3.04, 2.18, 8.17, 2.63, 7.69, 3.07, 3.7, 2.58,
2.09, 2.18, 2.16, 1.7, 6.8, 1.95, 4.28, 2.58, 2.09, 2.98,
2.16, 1.7, 6.8, 1.95, 4.28, 2.66, 3.14, 2.18, 8.17, 2.63,
7.69, 3.07, 3.7, 2.66, 3.14, 2.92, 2.16, 1.7, 6.8,
1.95, 4.28),
pH = c(9.94, 10.06, 9.98, 9.07, 8.99, 8.93, 9.24, 9.32, 9.94,
10, 9.98, 9.07, 8.99, 8.93, 9.24, 9.32, 10.04, 10, 10,
9.06, 8.99, 8.92, 9.18, 9.32, 10.04, 9.96, 9.98, 9.06,
8.99, 8.92, 9.18, 9.32, 9.94, 9.96, 9.98, 9.07, 8.99,
8.93, 9.24, 9.32, 9.94, 9.95, 10.09, 9.07, 8.99, 8.92,
9.18, 9.32, 10.04, 9.95, 10.09, 9.04, 9, 9.53, 9.31, 9.14,
10.04, 9.96, 9.97, 9.04, 9, 9.53, 9.31, 9.14, 10.06,
9.95, 9.97, 8.96, 9.02, 9.51, 9.39, 8.85, 10.06, 9.95,
10.09, 8.96, 9.02, 9.51, 9.39, 8.85, 10, 10, 9.97, 9.04, 9,
9.53, 9.31, 9.14, 10, 10, 9.06, 8.96, 9.02, 9.51, 9.39,
8.85),
Temperature = c(24.48, 24.95, 24.01, 26.67, 24.01, 23.5, 24.46, 24.65,
24.48, 24.92, 24.01, 26.67, 24.01, 23.5, 24.46, 24.65,
24.35, 24.92, 24, 26.46, 23.99, 23.43, 24.14, 24.67,
24.35, 24.04, 24.01, 26.46, 23.99, 23.43, 24.14, 24.67,
24.48, 24.04, 24.01, 26.67, 24.01, 23.5, 24.46, 24.65,
24.48, 24.01, 24.36, 26.67, 23.99, 23.43, 24.14, 24.67,
24.35, 24.01, 24.36, 24.88, 24.7, 26.14, 24.19, 22.65,
24.35, 24.04, 24.4, 24.88, 24.7, 26.14, 24.19, 22.65, 24.95,
24.01, 24.4, 24.3, 24.72, 26.15, 23.9, 20.65, 24.95,
24.01, 24.36, 24.3, 24.72, 26.15, 23.9, 20.65, 24.92, 24,
24.4, 24.88, 24.7, 26.14, 24.19, 22.65, 24.92, 24,
26.46, 24.3, 24.72, 26.15, 23.9, 20.65),
Turbidity = c(87.3, 95.2, 98.8, 24.9, 43.1, 49.1, 68.6, 69.8, 87.3,
113, 98.8, 24.9, 43.1, 49.1, 68.6, 69.8, 87.7, 113,
96.9, 39.7, 49, 47.6, 66.6, 74.1, 87.7, 98.5, 98.8, 39.7,
49, 47.6, 66.6, 74.1, 87.3, 98.5, 98.8, 24.9, 43.1,
49.1, 68.6, 69.8, 87.3, 95.9, 97.2, 24.9, 49, 47.6, 66.6,
74.1, 87.7, 95.9, 97.2, 46.9, 37.2, 74.2, 69.9, 71.3,
87.7, 98.5, 100.6, 46.9, 37.2, 74.2, 69.9, 71.3, 95.2, 95.9,
100.6, 42, 36.1, 74.4, 44.4, 71.3, 95.2, 95.9, 97.2,
42, 36.1, 74.4, 44.4, 71.3, 113, 96.9, 100.6, 46.9, 37.2,
74.2, 69.9, 71.3, 113, 96.9, 39.7, 42, 36.1, 74.4, 44.4,
71.3),
ORP = c(97, 108.4, 102.9, 94, 63.9, 64.7, 40.3, 48.1, 97,
111.9, 102.9, 94, 63.9, 64.7, 40.3, 48.1, 91.3, 111.9,
108, 23.7, 64.1, 60.4, 40.6, 49.8, 91.3, 98.9, 102.9,
23.7, 64.1, 60.4, 40.6, 49.8, 97, 98.9, 102.9, 94, 63.9,
64.7, 40.3, 48.1, 97, 92.8, 106.9, 94, 64.1, 60.4, 40.6,
49.8, 91.3, 92.8, 106.9, 90.4, 74.9, 67.7, 42, 37.1,
91.3, 98.9, 102.7, 90.4, 74.9, 67.7, 42, 37.1, 108.4,
92.8, 102.7, 53, 63.4, 73.9, 42, 37.5, 108.4, 92.8, 106.9,
53, 63.4, 73.9, 42, 37.5, 111.9, 108, 102.7, 90.4, 74.9,
67.7, 42, 37.1, 111.9, 111.9, 23.7, 53, 63.4, 73.9, 42,
37.5),
Ammonium = c(4.46, 4.62, 4.37, 5.6, 3.38, 3.21, 1.93, 2.21, 4.46,
4.75, 4.37, 5.6, 3.38, 3.21, 1.93, 2.21, 4.5, 4.75,
4.48, 7.64, 3.18, 3.09, 2.19, 2.14, 4.5, 3.97, 4.37, 7.64,
3.18, 3.09, 2.19, 2.14, 4.46, 3.97, 4.37, 5.6, 3.38,
3.21, 1.93, 2.21, 4.46, 4.06, 4.43, 5.6, 3.18, 3.09, 2.19,
2.14, 4.5, 4.06, 4.43, 2.81, 4.2, 4.5, 1.99, 2.07, 4.5,
3.97, 5.47, 2.81, 4.2, 4.5, 1.99, 2.07, 4.62, 4.06,
5.47, 3.19, 3.97, 3.3, 3.38, 1.18, 4.62, 4.06, 4.43, 3.19,
3.97, 3.3, 3.38, 1.18, 4.75, 4.48, 5.47, 2.81, 4.2, 4.5,
1.99, 2.07, 4.75, 4.48, 7.64, 3.19, 3.97, 3.3, 3.38,
1.18),
Nitrates = c(3.6, 3.98, 3.57, 3.53, 2.55, 2.28, 0.71, 1.2, 3.6,
3.21, 3.57, 3.53, 2.55, 2.28, 0.71, 1.2, 3.36, 3.21,
3.76, 3.66, 2.44, 2.06, 0.67, 0.87, 3.36, 4.35, 3.57,
3.66, 2.44, 2.06, 0.67, 0.87, 3.6, 4.35, 3.57, 3.53, 2.55,
2.28, 0.71, 1.2, 3.6, 3.68, 3.34, 3.53, 2.44, 2.06,
0.67, 0.87, 3.36, 3.68, 3.34, 2.25, 2.99, 1.3, 0.74, 0.68,
3.36, 4.35, 4.11, 2.25, 2.99, 1.3, 0.74, 0.68, 3.98,
3.68, 4.11, 2.1, 2.55, 1.11, 0.69, 0.69, 3.98, 3.68, 3.34,
2.1, 2.55, 1.11, 0.69, 0.69, 3.21, 3.76, 4.11, 2.25,
2.99, 1.3, 0.74, 0.68, 3.21, 3.76, 3.66, 2.1, 2.55, 1.11,
0.69, 0.69),
BGA = c(279133, 279113, 265304, 11603, 218471, 226658, 267389,
267144, 279133, 270193, 265304, 11603, 218471, 226658,
267389, 267144, 267716, 270193, 279119, 66608, 213497,
242331, 267007, 261363, 267716, 279109, 265304, 66608,
213497, 242331, 267007, 261363, 279133, 279109, 265304,
11603, 218471, 226658, 267389, 267144, 279133, 279100,
279099, 11603, 213497, 242331, 267007, 261363, 267716,
279100, 279099, 246979, 169175, 272332, 269519, 272096,
267716, 279109, 278137, 246979, 169175, 272332, 269519,
272096, 279113, 279100, 278137, 174271, 159046, 271503,
261124, 271303, 279113, 279100, 279099, 174271, 159046,
271503, 261124, 271303, 270193, 279119, 278137, 246979,
169175, 272332, 269519, 272096, 270193, 279119, 66608,
174271, 159046, 271503, 261124, 271303),
Chlrophyll = c(36.6, 38.9, 40, 44, 35, 32.8, 44, 44.8, 36.6, 39.8,
40, 44, 35, 32.8, 44, 44.8, 39.6, 39.8, 38.5, 45.1,
42.2, 50.2, 43.6, 45.1, 39.6, 43, 40, 45.1, 42.2, 50.2,
43.6, 45.1, 36.6, 43, 40, 44, 35, 32.8, 44, 44.8, 36.6,
50, 43, 44, 42.2, 50.2, 43.6, 45.1, 39.6, 50, 43, 34,
38.4, 58.5, 44.3, 40, 39.6, 43, 53.8, 34, 38.4, 58.5, 44.3,
40, 38.9, 50, 53.8, 50.3, 55.9, 51.3, 68.2, 34.2, 38.9,
50, 43, 50.3, 55.9, 51.3, 68.2, 34.2, 39.8, 38.5, 53.8,
34, 38.4, 58.5, 44.3, 40, 39.8, 38.5, 45.1, 50.3, 55.9,
51.3, 68.2, 34.2)
)
library(readxl)
Genera <- read_excel("~/RSTUDIO/Datos_cianobacterias/Ciano_genera.xlsx")
#NORMALIZATION OF RAW READ TABLE LOG2
#OTUS en filas, columnas son las muestras, es al reves que para el resto de pasos
data<- Genera
replicates <- as.data.frame(colnames(data)[-1])
colnames(replicates) <- "replicates"
attach(Genera)
rwnames <- OTU
data <- as.matrix(data[,-1])
rownames(data) <- rwnames
data[is.na(data)] <- 0
library(DESeq2)
dds <- DESeqDataSetFromMatrix(data, replicates, ~ replicates)
cts <- counts(dds)
geoMeans <- apply(cts, 1, function(row) if (all(row == 0)) 0 else exp(mean(log(row[row != 0]))))
dds <- estimateSizeFactors(dds, geoMeans=geoMeans)
norm <- counts(dds, normalized=TRUE)
logcounts <- log2( counts(dds, normalized=TRUE) + 1 )# in case you want a log transformed normalized count
write.csv(logcounts, file = "~/RSTUDIO/Datos_cianobacterias/ciano_genero_normalizados.csv")
#log counts are like variance stabilized counts (https://support.bioconductor.org/p/76712/)
#DATA ANALYSIS
#Introduce two tables, one has to have the samples in the first column and the
#species in the first line with the reads per sample.
#The second table is with environmental parameters, it has to have the samples
#in the first column and the environmental parameters in the first line.
#You can have the environmental parameters as presence/absence with 1 and 0 or
#with numbers.
library(readxl)
Genera <- read_excel("~/RSTUDIO/Datos_cianobacterias/ciano_genera_normalizados_traspu.xlsx")
dim(Genera)
library(vegan)
# Solo acepta variables numéricas así que debes obviar las variables de texto
rda_object <- vegan::rda(Genera[,-1], metadata[,7:15], scale=TRUE)
plot(rda_object)
library(BiodiversityR)
PCAsignificance(rda_object, axes=8)
library(devtools)
library(ggord)
levels(metadata$Month) <- c("July","August","September")
Genera <- metadata$Month
bg <- c("#ff7f00","#1f78b4","#ffff33","#a6cee3","#33a02c","#e31a1c") # 6 nice colors for our ecotypes
# axes 1 & 2
plot(rda_object, type = "n", scaling=3, display = "sites", xlab = "RDA1 (14.59%)", ylab = "RDA2 (9.68%)")
points(rda_object, display="species", pch=20, cex=0.7, col="gray32", scaling=3) # the SNPs
points(rda_object, display="sites", pch=21, cex=1.3, col="gray32", scaling=3, bg=bg[Genera]) # the wolves
text(rda_object, scaling=3, display="bp", col="#0868ac", cex=1) # the predictors
legend("bottomright", legend=levels(Genera), bty="n", col="gray32", pch=21, cex=1, pt.bg=bg)