MicrobiotaProcess Group2 vs Group6 (v1)

gene_x 0 like s 1024 view s

Tags: plot, packages, R, processing, pipeline

https://bioconductor.org/packages/release/bioc/vignettes/MicrobiotaProcess/inst/doc//MicrobiotaProcess.html

diff_analysis_Group2_vs_Group6

PCoA2_Group2_vs_Group6

PCoA2_Group2_vs_Group6

  1. code of MicrobiotaProcess_Group2_vs_Group6.R

    1. # https://bioconductor.org/packages/release/bioc/vignettes/MicrobiotaProcess/inst/doc//MicrobiotaProcess.html
    2. # -----------------------------------
    3. # ---- prepare the R environment ----
    4. #Rscript MicrobiotaProcess.R
    5. #NOTE: exit R script, then login again R-environment; rm -rf Phyloseq*_cache
    6. rmarkdown::render('Phyloseq.Rmd',output_file='Phyloseq.html')
    7. # -----------------------------
    8. # ---- 3.1. bridges other tools
    9. ##https://github.com/YuLab-SMU/MicrobiotaProcess
    10. ##https://www.bioconductor.org/packages/release/bioc/vignettes/MicrobiotaProcess/inst/doc/MicrobiotaProcess.html
    11. ##https://chiliubio.github.io/microeco_tutorial/intro.html#framework
    12. ##https://yiluheihei.github.io/microbiomeMarker/reference/plot_cladogram.html
    13. #BiocManager::install("MicrobiotaProcess")
    14. #install.packages("microeco")
    15. #install.packages("ggalluvial")
    16. #install.packages("ggh4x")
    17. library(MicrobiotaProcess)
    18. library(microeco)
    19. library(ggalluvial)
    20. library(ggh4x)
    21. library(gghalves)
    22. ## Convert the phyloseq object to a MicrobiotaProcess object
    23. #mp <- as.MicrobiotaProcess(ps.ng.tax)
    24. #mt <- phyloseq2microeco(ps.ng.tax) #--> ERROR
    25. #abundance_table <- mt$abun_table
    26. #taxonomy_table <- mt$tax_table
    27. #ps.ng.tax_abund <- phyloseq::filter_taxa(ps.ng.tax, function(x) sum(x > total*0.01) > 0, TRUE)
    28. #ps.ng.tax_most = phyloseq::filter_taxa(ps.ng.tax_rel, function(x) mean(x) > 0.001, TRUE)
    29. ##OPTION1 (NOT_USED): take all samples, prepare ps.ng.tax_abund --> mpse_abund
    30. ##mpse <- ps.ng.tax %>% as.MPSE()
    31. #mpse_abund <- ps.ng.tax_abund %>% as.MPSE()
    32. ##OPTION2 (USED!): take partial samples, prepare ps.ng.tax or ps.ng.tax_abund (2 replacements!)--> ps.ng.tax_sel --> mpse_abund
    33. ps.ng.tax_sel <- ps.ng.tax_abund
    34. ##otu_table(ps.ng.tax_sel) <- otu_table(ps.ng.tax)[,c("1","2","5","6","7", "15","16","17","18","19","20", "29","30","31","32", "40","41","42","43","44","46")]
    35. ##NOTE: Only choose Group2, Group4, Group6, Group8
    36. #> ps.ng.tax_sel
    37. #otu_table() OTU Table: [ 37465 taxa and 29 samples ]
    38. #sample_data() Sample Data: [ 29 samples by 10 sample variables ]
    39. #tax_table() Taxonomy Table: [ 37465 taxa by 7 taxonomic ranks ]
    40. #phy_tree() Phylogenetic Tree: [ 37465 tips and 37461 internal nodes ]
    41. #-Group4: "21","22","23","24","25","26","27","28",
    42. #-Group8: , "47","48","49","50","52","53","55"
    43. otu_table(ps.ng.tax_sel) <- otu_table(ps.ng.tax_abund)[,c("8","9","10","12","13","14", "33","34","35","36","37","38","39","51")]
    44. mpse_abund <- ps.ng.tax_sel %>% as.MPSE()
    45. # A MPSE-tibble (MPSE object) abstraction: 2,352 × 20
    46. # NOTE mpse_abund contains 20 variables: OTU, Sample, Abundance, BarcodeSequence, LinkerPrimerSequence, FileInput, Group,
    47. # Sex_age <chr>, pre_post_stroke <chr>, Conc <dbl>, Vol_50ng <dbl>, Vol_PCR <dbl>, Description <chr>,
    48. # Domain <chr>, Phylum <chr>, Class <chr>, Order <chr>, Family <chr>, Genus <chr>, Species <chr>
    49. # -----------------------------------
    50. # ---- 3.2. alpha diversity analysis
    51. # Rarefied species richness + RareAbundance
    52. mpse_abund %<>% mp_rrarefy()
    53. # 'chunks' represent the split number of each sample to calculate alpha
    54. # diversity, default is 400. e.g. If a sample has total 40000
    55. # reads, if chunks is 400, it will be split to 100 sub-samples
    56. # (100, 200, 300,..., 40000), then alpha diversity index was
    57. # calculated based on the sub-samples.
    58. # '.abundance' the column name of abundance, if the '.abundance' is not be
    59. # rarefied calculate rarecurve, user can specific 'force=TRUE'.
    60. mpse_abund %<>%
    61. mp_cal_rarecurve(
    62. .abundance = RareAbundance,
    63. chunks = 400
    64. )
    65. # The RareAbundanceRarecurve column will be added the colData slot
    66. # automatically (default action="add")
    67. #NOTE mpse_abund contains 22 varibles = 20 varibles + RareAbundance <dbl> + RareAbundanceRarecurve <list>
    68. # default will display the confidence interval around smooth.
    69. # se=TRUE
    70. # NOTE that two colors #c("#00A087FF", "#3C5488FF") for .group = pre_post_stroke; four colors c("#1f78b4", "#33a02c", "#e31a1c", "#6a3d9a") for .group = Group;
    71. p1 <- mpse_abund %>%
    72. mp_plot_rarecurve(
    73. .rare = RareAbundanceRarecurve,
    74. .alpha = Observe,
    75. )
    76. p2 <- mpse_abund %>%
    77. mp_plot_rarecurve(
    78. .rare = RareAbundanceRarecurve,
    79. .alpha = Observe,
    80. .group = Group
    81. ) +
    82. scale_color_manual(values=c("#1f78b4", "#e31a1c")) +
    83. scale_fill_manual(values=c("#1f78b4", "#e31a1c"), guide="none")
    84. # combine the samples belong to the same groups if plot.group=TRUE
    85. p3 <- mpse_abund %>%
    86. mp_plot_rarecurve(
    87. .rare = RareAbundanceRarecurve,
    88. .alpha = "Observe",
    89. .group = Group,
    90. plot.group = TRUE
    91. ) +
    92. scale_color_manual(values=c("#1f78b4", "#e31a1c")) +
    93. scale_fill_manual(values=c("#1f78b4", "#e31a1c"),guide="none")
    94. png("rarefaction_of_samples_or_groups.png", width=1080, height=600)
    95. p1 + p2 + p3
    96. dev.off()
    97. # ------------------------------------------
    98. # 3.3. calculate alpha index and visualization
    99. library(ggplot2)
    100. library(MicrobiotaProcess)
    101. mpse_abund %<>%
    102. mp_cal_alpha(.abundance=RareAbundance)
    103. mpse_abund
    104. #NOTE mpse_abund contains 28 varibles = 22 varibles + Observe <dbl>, Chao1 <dbl>, ACE <dbl>, Shannon <dbl>, Simpson <dbl>, Pielou <dbl>
    105. f1 <- mpse_abund %>%
    106. mp_plot_alpha(
    107. .group=Group,
    108. .alpha=c(Observe, Chao1, ACE, Shannon, Simpson, Pielou)
    109. ) +
    110. scale_fill_manual(values=c("#1f78b4", "#e31a1c"), guide="none") +
    111. scale_color_manual(values=c("#1f78b4", "#e31a1c"), guide="none")
    112. f2 <- mpse_abund %>%
    113. mp_plot_alpha(
    114. .alpha=c(Observe, Chao1, ACE, Shannon, Simpson, Pielou)
    115. )
    116. #ps.ng.tax_sel contais only pre samples --> f1 cannot be generated!
    117. png("alpha_diversity_comparison.png", width=1400, height=600)
    118. f1 / f2
    119. dev.off()
    120. # -------------------------------------------
    121. # 3.4. The visualization of taxonomy abundance (Class)
    122. mpse_abund %<>%
    123. mp_cal_abundance( # for each samples
    124. .abundance = RareAbundance
    125. ) %>%
    126. mp_cal_abundance( # for each groups
    127. .abundance=RareAbundance,
    128. .group=Group
    129. )
    130. mpse_abund
    131. #NOTE mpse_abund contains 29 varibles = 28 varibles + RelRareAbundanceBySample <dbl>
    132. # visualize the relative abundance of top 20 phyla for each sample.
    133. # .group=time,
    134. p1 <- mpse_abund %>%
    135. mp_plot_abundance(
    136. .abundance=RareAbundance,
    137. taxa.class = Class,
    138. topn = 20,
    139. relative = TRUE
    140. )
    141. # visualize the abundance (rarefied) of top 20 phyla for each sample.
    142. # .group=time,
    143. p2 <- mpse_abund %>%
    144. mp_plot_abundance(
    145. .abundance=RareAbundance,
    146. taxa.class = Class,
    147. topn = 20,
    148. relative = FALSE
    149. )
    150. png("relative_abundance_and_abundance.png", width= 1200, height=600) #NOT PRODUCED!
    151. p1 / p2
    152. dev.off()
    153. #----
    154. h1 <- mpse_abund %>%
    155. mp_plot_abundance(
    156. .abundance = RareAbundance,
    157. .group = Group,
    158. taxa.class = Class,
    159. relative = TRUE,
    160. topn = 20,
    161. geom = 'heatmap',
    162. features.dist = 'euclidean',
    163. features.hclust = 'average',
    164. sample.dist = 'bray',
    165. sample.hclust = 'average'
    166. )
    167. h2 <- mpse_abund %>%
    168. mp_plot_abundance(
    169. .abundance = RareAbundance,
    170. .group = Group,
    171. taxa.class = Class,
    172. relative = FALSE,
    173. topn = 20,
    174. geom = 'heatmap',
    175. features.dist = 'euclidean',
    176. features.hclust = 'average',
    177. sample.dist = 'bray',
    178. sample.hclust = 'average'
    179. )
    180. # the character (scale or theme) of figure can be adjusted by set_scale_theme
    181. # refer to the mp_plot_dist
    182. png("relative_abundance_and_abundance_heatmap.png", width= 1200, height=600)
    183. aplot::plot_list(gglist=list(h1, h2), tag_levels="A")
    184. dev.off()
    185. # visualize the relative abundance of top 20 class for each .group (Group)
    186. p3 <- mpse_abund %>%
    187. mp_plot_abundance(
    188. .abundance=RareAbundance,
    189. .group=Group,
    190. taxa.class = Class,
    191. topn = 20,
    192. plot.group = TRUE
    193. )
    194. # visualize the abundance of top 20 phyla for each .group (time)
    195. p4 <- mpse_abund %>%
    196. mp_plot_abundance(
    197. .abundance=RareAbundance,
    198. .group= Group,
    199. taxa.class = Class,
    200. topn = 20,
    201. relative = FALSE,
    202. plot.group = TRUE
    203. )
    204. png("relative_abundance_and_abundance_groups.png", width= 1000, height=1000)
    205. p3 / p4
    206. dev.off()
    207. # ---------------------------
    208. # 3.5. Beta diversity analysis
    209. # ---------------------------------------------
    210. # 3.5.1 The distance between samples or groups
    211. # standardization
    212. # mp_decostand wraps the decostand of vegan, which provides
    213. # many standardization methods for community ecology.
    214. # default is hellinger, then the abundance processed will
    215. # be stored to the assays slot.
    216. mpse_abund %<>%
    217. mp_decostand(.abundance=Abundance)
    218. mpse_abund
    219. #NOTE mpse_abund contains 30 varibles = 29 varibles + hellinger <dbl>
    220. # calculate the distance between the samples.
    221. # the distance will be generated a nested tibble and added to the
    222. # colData slot.
    223. mpse_abund %<>% mp_cal_dist(.abundance=hellinger, distmethod="bray")
    224. mpse_abund
    225. #NOTE mpse_abund contains 31 varibles = 30 varibles + bray <list>
    226. # mp_plot_dist provides there methods to visualize the distance between the samples or groups
    227. # when .group is not provided, the dot heatmap plot will be return
    228. p1 <- mpse_abund %>% mp_plot_dist(.distmethod = bray)
    229. png("distance_between_samples.png", width= 1000, height=1000)
    230. p1
    231. dev.off()
    232. # when .group is provided, the dot heatmap plot with group information will be return.
    233. p2 <- mpse_abund %>% mp_plot_dist(.distmethod = bray, .group = Group)
    234. # The scale or theme of dot heatmap plot can be adjusted using set_scale_theme function.
    235. p2 %>% set_scale_theme(
    236. x = scale_fill_manual(
    237. values=c("#1f78b4", "#e31a1c"), #c("orange", "deepskyblue"),
    238. guide = guide_legend(
    239. keywidth = 1,
    240. keyheight = 0.5,
    241. title.theme = element_text(size=8),
    242. label.theme = element_text(size=6)
    243. )
    244. ),
    245. aes_var = Group # specific the name of variable
    246. ) %>%
    247. set_scale_theme(
    248. x = scale_color_gradient(
    249. guide = guide_legend(keywidth = 0.5, keyheight = 0.5)
    250. ),
    251. aes_var = bray
    252. ) %>%
    253. set_scale_theme(
    254. x = scale_size_continuous(
    255. range = c(0.1, 3),
    256. guide = guide_legend(keywidth = 0.5, keyheight = 0.5)
    257. ),
    258. aes_var = bray
    259. )
    260. png("distance_between_samples_with_group_info.png", width= 1000, height=1000)
    261. p2
    262. dev.off()
    263. # when .group is provided and group.test is TRUE, the comparison of different groups will be returned
    264. # Assuming p3 is a ggplot object after mp_plot_dist call
    265. p3 <- mpse_abund %>%
    266. mp_plot_dist(.distmethod = bray, .group = Group, group.test = TRUE, textsize = 6) +
    267. theme(
    268. axis.title.x = element_text(size = 14), # Customize x-axis label face = "bold"
    269. axis.title.y = element_text(size = 14), # Customize y-axis label
    270. axis.text.x = element_text(size = 14), # Customize x-axis ticks
    271. axis.text.y = element_text(size = 14) # Customize y-axis ticks
    272. )
    273. # Save the plot with the new theme settings
    274. png("Comparison_of_Bray_Distances.png", width = 1000, height = 1000)
    275. print(p3) # Ensure that p3 is explicitly printed in the device
    276. dev.off()
    277. # Extract Bray-Curtis Distance Values and save them in a Excel-table.
    278. library(dplyr)
    279. library(tidyr)
    280. library(openxlsx)
    281. # Define the sample numbers vector
    282. sample_numbers <- c("8", "9", "10", "12", "13", "14", "33", "34", "35", "36", "37", "38", "39", "51")
    283. # Consolidate the list of tibbles using the actual sample numbers
    284. bray_data <- bind_rows(
    285. lapply(seq_along(mpse_abund$bray), function(i) {
    286. tibble(
    287. Sample1 = sample_numbers[i], # Use actual sample number
    288. Sample2 = mpse_abund$bray[[i]]$braySampley,
    289. BrayDistance = mpse_abund$bray[[i]]$bray
    290. )
    291. }),
    292. .id = "PairID"
    293. )
    294. # Print the data frame to check the output
    295. print(bray_data)
    296. # Write the data frame to an Excel file
    297. write.xlsx(bray_data, file = "Bray_Curtis_Distances.xlsx")
    298. #DELETE the column "PairID" in Excel file
    299. # -----------------------
    300. # 3.5.2 The PCoA analysis
    301. #install.packages("corrr")
    302. library(corrr)
    303. #install.packages("ggside")
    304. library(ggside)
    305. mpse_abund %<>%
    306. mp_cal_pcoa(.abundance=hellinger, distmethod="bray")
    307. # The dimensions of ordination analysis will be added the colData slot (default).
    308. mpse_abund
    309. mpse_abund %>% print(width=380, n=2)
    310. #NOTE mpse_abund contains 34 varibles = 31 varibles + `PCo1 (30.16%)` <dbl>, `PCo2 (15.75%)` <dbl>, `PCo3 (10.53%)` <dbl>
    311. #BUG why 36 variables in mpse_abund %>% print(width=380, n=1) [RareAbundanceBySample <list>, RareAbundanceByGroup <list>]
    312. #> methods(class=class(mpse_abund))
    313. # [1] [ [[<- [<-
    314. # [4] $ $<- arrange
    315. # [7] as_tibble as.data.frame as.phyloseq
    316. #[10] coerce coerce<- colData<-
    317. #[13] distinct filter group_by
    318. #[16] left_join mp_adonis mp_aggregate_clade
    319. #[19] mp_aggregate mp_anosim mp_balance_clade
    320. #[22] mp_cal_abundance mp_cal_alpha mp_cal_cca
    321. #[25] mp_cal_clust mp_cal_dca mp_cal_dist
    322. #[28] mp_cal_nmds mp_cal_pca mp_cal_pcoa
    323. #[31] mp_cal_pd_metric mp_cal_rarecurve mp_cal_rda
    324. #[34] mp_cal_upset mp_cal_venn mp_decostand
    325. #[37] mp_diff_analysis mp_diff_clade mp_envfit
    326. #[40] mp_extract_abundance mp_extract_assays mp_extract_dist
    327. #[43] mp_extract_feature mp_extract_internal_attr mp_extract_rarecurve
    328. #[46] mp_extract_refseq mp_extract_sample mp_extract_taxonomy
    329. #[49] mp_extract_tree mp_filter_taxa mp_mantel
    330. #[52] mp_mrpp mp_plot_abundance mp_plot_alpha
    331. #[55] mp_plot_diff_boxplot mp_plot_diff_res mp_plot_dist
    332. #[58] mp_plot_ord mp_plot_rarecurve mp_plot_upset
    333. #[61] mp_plot_venn mp_rrarefy mp_select_as_tip
    334. #[64] mp_stat_taxa mutate otutree
    335. #[67] otutree<- print pull
    336. #[70] refsequence refsequence<- rename
    337. #[73] rownames<- select show
    338. # [ reached getOption("max.print") -- omitted 6 entries ]
    339. #see '?methods' for accessing help and source code
    340. # We also can perform adonis or anosim to check whether it is significant to the dissimilarities of groups.
    341. mpse_abund %<>%
    342. mp_adonis(.abundance=hellinger, .formula=~Group, distmethod="bray", permutations=9999, action="add")
    343. mpse_abund %>% mp_extract_internal_attr(name=adonis)
    344. #NOTE mpse_abund contains 34 varibles, no new variable, it has been saved in mpse_abund and can be extracted with "mpse_abund %>% mp_extract_internal_attr(name='adonis')"
    345. #The result of adonis has been saved to the internal attribute !
    346. #It can be extracted using this-object %>% mp_extract_internal_attr(name='adonis')
    347. #The object contained internal attribute: PCoA ADONIS
    348. #Permutation test for adonis under reduced model
    349. #Terms added sequentially (first to last)
    350. #Permutation: free
    351. #Number of permutations: 9999
    352. #
    353. #vegan::adonis2(formula = .formula, data = sampleda, permutations = permutations, method = distmethod)
    354. # Df SumOfSqs R2 F Pr(>F)
    355. #Group 1 0.23448 0.22659 3.5158 5e-04 ***
    356. #Residual 12 0.80032 0.77341
    357. #Total 13 1.03480 1.00000
    358. #---
    359. #Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
    360. # ("1","2","5","6","7", "15","16","17","18","19","20", "29","30","31","32", "40","41","42","43","44","46")
    361. #div.df2[div.df2 == "Group1"] <- "aged.post"
    362. #div.df2[div.df2 == "Group3"] <- "young.post"
    363. #div.df2[div.df2 == "Group5"] <- "aged.post"
    364. #div.df2[div.df2 == "Group7"] <- "young.post"
    365. # ("8","9","10","12","13","14", "21","22","23","24","25","26","27","28", "33","34","35","36","37","38","39","51", "47","48","49","50","52","53","55")
    366. #div.df2[div.df2 == "Group2"] <- "aged.pre"
    367. #div.df2[div.df2 == "Group4"] <- "young.pre"
    368. #div.df2[div.df2 == "Group6"] <- "aged.pre"
    369. #div.df2[div.df2 == "Group8"] <- "young.pre"
    370. #Group1: f.aged and post
    371. #Group2: f.aged and pre
    372. #Group3: f.young and post
    373. #Group4: f.young and pre
    374. #Group5: m.aged and post
    375. #Group6: m.aged and pre
    376. #Group7: m.young and post
    377. #Group8: m.young and pre
    378. #[,c("1","2","5","6","7", "8","9","10","12","13","14")]
    379. #[,c("15","16","17","18","19","20", "21","22","23","24","25","26","27","28")]
    380. #[,c("29","30","31","32", "33","34","35","36","37","38","39","51")]
    381. #[,c("40","41","42","43","44","46", "47","48","49","50","52","53","55")]
    382. #For the first set:
    383. #a6cee3: This is a light blue color, somewhat pastel and soft.
    384. #b2df8a: A soft, pale green, similar to a light lime.
    385. #fb9a99: A soft pink, slightly peachy or salmon-like.
    386. #cab2d6: A pale purple, reminiscent of lavender or a light mauve.
    387. #For the second set:
    388. #1f78b4: This is a strong, vivid blue, close to cobalt or a medium-dark blue.
    389. #33a02c: A medium forest green, vibrant and leafy.
    390. #e31a1c: A bright red, very vivid, similar to fire engine red.
    391. #6a3d9a: This would be described as a deep purple, akin to a dark lavender or plum.
    392. p1 <- mpse_abund %>%
    393. mp_plot_ord(
    394. .ord = pcoa,
    395. .group = Group,
    396. .color = Group,
    397. .size = 2.4,
    398. .alpha = 1,
    399. ellipse = TRUE,
    400. show.legend = FALSE # don't display the legend of stat_ellipse
    401. ) +
    402. scale_fill_manual(
    403. #values = c("#a6cee3", "#1f78b4", "#b2df8a", "#33a02c", "#fb9a99", "#e31a1c", "#cab2d6", "#6a3d9a"),
    404. #values = c("#a6cee3", "#b2df8a", "#fb9a99", "#cab2d6"),
    405. values = c("#1f78b4", "#e31a1c"),
    406. guide = guide_legend(keywidth=1.6, keyheight=1.6, label.theme=element_text(size=12))
    407. ) +
    408. scale_color_manual(
    409. #values=c("#a6cee3", "#1f78b4", "#b2df8a", "#33a02c", "#fb9a99", "#e31a1c", "#cab2d6", "#6a3d9a"),
    410. #values = c("#a6cee3", "#b2df8a", "#fb9a99", "#cab2d6"),
    411. values = c("#1f78b4", "#e31a1c"),
    412. guide = guide_legend(keywidth=1.6, keyheight=1.6, label.theme=element_text(size=12))
    413. )
    414. #scale_fill_manual(values=c("#FF0000", "#000000", "#0000FF", "#C0C0C0", "#00FF00", "#FFFF00", "#00FFFF", "#FFA500")) +
    415. #scale_color_manual(values=c("#FF0000", "#000000", "#0000FF", "#C0C0C0", "#00FF00", "#FFFF00", "#00FFFF", "#FFA500"))
    416. #scale_fill_manual(values=c("#00A087FF", "#3C5488FF")) +
    417. #scale_color_manual(values=c("#00A087FF", "#3C5488FF"))
    418. #png("PCoA.png", width= 1000, height=1000)
    419. #svg("PCoA.svg", width= 11, height=10)
    420. #svg("PCoA_.svg", width=10, height=10)
    421. #svg("PCoA.svg")
    422. pdf("PCoA_Group2_vs_Group6.pdf")
    423. p1
    424. dev.off()
    425. #FF0000: Red
    426. #000000: Black
    427. #0000FF: Blue
    428. #C0C0C0: Silver
    429. #00FF00: Lime (often referred to simply as Green in web colors)
    430. #FFFF00: Yellow
    431. #00FFFF: Aqua (also known as Cyan)
    432. #FFA500: Orange
    433. # The size of point also can be mapped to other variables such as Observe, or Shannon
    434. # Then the alpha diversity and beta diversity will be displayed simultaneously.
    435. p2 <- mpse_abund %>%
    436. mp_plot_ord(
    437. .ord = pcoa,
    438. .group = Group,
    439. .color = Group,
    440. .size = Shannon,
    441. .alpha = Observe,
    442. ellipse = TRUE,
    443. show.legend = FALSE # don't display the legend of stat_ellipse
    444. ) +
    445. scale_fill_manual(
    446. values = c("#1f78b4", "#e31a1c"), #only needs four colors.
    447. #values = c("#FF0000", "#000000", "#0000FF", "#C0C0C0", "#00FF00", "#FFFF00", "#00FFFF", "#FFA500"),
    448. guide = guide_legend(keywidth=0.6, keyheight=0.6, label.theme=element_text(size=8))
    449. ) +
    450. scale_color_manual(
    451. values = c("#1f78b4", "#e31a1c"), #only needs four colors.
    452. #values=c("#FF0000", "#000000", "#0000FF", "#C0C0C0", "#00FF00", "#FFFF00", "#00FFFF", "#FFA500"),
    453. guide = guide_legend(keywidth=0.6, keyheight=0.6, label.theme=element_text(size=8))
    454. ) +
    455. scale_size_continuous(
    456. range=c(0.5, 3),
    457. guide = guide_legend(keywidth=0.6, keyheight=0.6, label.theme=element_text(size=8))
    458. )
    459. pdf("PCoA2_Group2_vs_Group6.pdf")
    460. p2
    461. dev.off()
    462. # ------------------------------------------
    463. # 3.5.3 Hierarchical cluster (tree) analysis
    464. #input should contain hellinger!
    465. mpse_abund %<>%
    466. mp_cal_clust(
    467. .abundance = hellinger,
    468. distmethod = "bray",
    469. hclustmethod = "average", # (UPGAE)
    470. action = "add" # action is used to control which result will be returned
    471. )
    472. mpse_abund
    473. mpse_abund %>% print(width=380, n=2)
    474. #NOTE mpse_abund contains 34 varibles, no new variable, the column bray has been new calculated!
    475. # if action = 'add', the result of hierarchical cluster will be added to the MPSE object
    476. # mp_extract_internal_attr can extract it. It is a treedata object, so it can be visualized
    477. # by ggtree.
    478. sample.clust <- mpse_abund %>% mp_extract_internal_attr(name='SampleClust')
    479. #The object contained internal attribute: PCoA ADONIS SampleClust
    480. sample.clust
    481. #--> The associated data tibble abstraction: 27 × 30
    482. library(ggtree)
    483. p <- ggtree(sample.clust) +
    484. geom_tippoint(aes(color=Group)) +
    485. geom_tiplab(as_ylab = TRUE) +
    486. ggplot2::scale_x_continuous(expand=c(0, 0.01))
    487. png("hierarchical_cluster1.png", width= 1000, height=800)
    488. p
    489. dev.off()
    490. #https://bioconductor.org/packages/release/bioc/vignettes/MicrobiotaProcess/inst/doc//MicrobiotaProcess.html
    491. # mapping = aes(x = RelRareAbundanceBySample-->Group,
    492. # y = Sample-->Group,
    493. # fill = Phyla
    494. # ),
    495. library(ggtreeExtra)
    496. library(ggplot2)
    497. # Extract relative abundance of phyla
    498. phyla.tb <- mpse_abund %>%
    499. mp_extract_abundance(taxa.class=Phylum, topn=30)
    500. # The abundance of each samples is nested, it can be flatted using the unnest of tidyr.
    501. phyla.tb %<>% tidyr::unnest(cols=RareAbundanceBySample) %>% dplyr::rename(Phyla="label")
    502. phyla.tb
    503. phyla.tb %>% print(width=380, n=10)
    504. p1 <- p +
    505. geom_fruit(
    506. data=phyla.tb,
    507. geom=geom_col,
    508. mapping = aes(x = RelRareAbundanceBySample,
    509. y = Sample,
    510. fill = Phyla
    511. ),
    512. orientation = "y",
    513. #offset = 0.4,
    514. pwidth = 3,
    515. axis.params = list(axis = "x",
    516. title = "The relative abundance of phyla (%)",
    517. title.size = 4,
    518. text.size = 2,
    519. vjust = 1),
    520. grid.params = list()
    521. )
    522. png("hierarchical_cluster2_Phyla.png", width = 1000, height = 800)
    523. p1
    524. dev.off()
    525. # Extract relative abundance of classes
    526. class.tb <- mpse_abund %>%
    527. mp_extract_abundance(taxa.class = Class, topn = 30)
    528. # Flatten and rename the columns
    529. class.tb %<>% tidyr::unnest(cols = RareAbundanceBySample) %>% dplyr::rename(Class = "label")
    530. # View the data frame
    531. class.tb
    532. # Create the plot
    533. p1 <- p +
    534. geom_fruit(
    535. data = class.tb,
    536. geom = geom_col,
    537. mapping = aes(x = RelRareAbundanceBySample,
    538. y = Sample,
    539. fill = Class
    540. ),
    541. orientation = "y",
    542. pwidth = 3,
    543. axis.params = list(axis = "x",
    544. title = "The relative abundance of classes (%)",
    545. title.size = 4,
    546. text.size = 2,
    547. vjust = 1),
    548. grid.params = list()
    549. )
    550. # Save the plot to a file #ERROR-->NEED to be DEBUGGED!
    551. png("hierarchical_cluster2_Class.png", width = 1000, height = 800)
    552. print(p1)
    553. dev.off()
    554. # -----------------------
    555. # 3.6 Biomarker discovery
    556. library(ggtree)
    557. library(ggtreeExtra)
    558. library(ggplot2)
    559. library(MicrobiotaProcess)
    560. library(tidytree)
    561. library(ggstar)
    562. library(forcats)
    563. mpse_abund %>% print(width=150)
    564. #mpse_abund %<>%
    565. # mp_cal_abundance( # for each samples
    566. # .abundance = RareAbundance
    567. # ) %>%
    568. # mp_cal_abundance( # for each groups
    569. # .abundance=RareAbundance,
    570. # .group=Group
    571. # )
    572. #mpse_abund
    573. mpse_abund %<>%
    574. mp_diff_analysis(
    575. .abundance = RelRareAbundanceBySample,
    576. .group = Group,
    577. first.test.alpha = 0.01, filter.p="pvalue"
    578. )
    579. # The result is stored to the taxatree or otutree slot, you can use mp_extract_tree to extract the specific slot.
    580. taxa.tree <- mpse_abund %>%
    581. mp_extract_tree(type="taxatree")
    582. taxa.tree
    583. ## And the result tibble of different analysis can also be extracted
    584. ## with tidytree (>=0.3.5)
    585. taxa.tree %>% select(label, nodeClass, LDAupper, LDAmean, LDAlower, Sign_Group, pvalue, fdr) %>% dplyr::filter(!is.na(fdr))
    586. taxa.tree %>% print(width=150, n=1)
    587. # Since taxa.tree is treedata object, it can be visualized by ggtree and ggtreeExtra
    588. p1 <- ggtree(
    589. taxa.tree,
    590. layout="radial",
    591. size = 0.3
    592. ) +
    593. geom_point(
    594. data = td_filter(!isTip),
    595. fill="white",
    596. size=1,
    597. shape=21
    598. )
    599. # display the high light of phylum clade.
    600. p2 <- p1 +
    601. geom_hilight(
    602. data = td_filter(nodeClass == "Phylum"),
    603. mapping = aes(node = node, fill = label)
    604. )
    605. # display the relative abundance of features(OTU)
    606. p3 <- p2 +
    607. ggnewscale::new_scale("fill") +
    608. geom_fruit(
    609. data = td_unnest(RareAbundanceBySample),
    610. geom = geom_star,
    611. mapping = aes(
    612. x = fct_reorder(Sample, Group, .fun=min),
    613. size = RelRareAbundanceBySample,
    614. fill = Group,
    615. subset = RelRareAbundanceBySample > 0
    616. ),
    617. starshape = 13,
    618. starstroke = 0.25,
    619. offset = 0.03,
    620. pwidth = 0.4,
    621. grid.params = list(linetype=2)
    622. ) +
    623. scale_size_continuous(
    624. name="Relative Abundance (%)",
    625. range = c(.5, 3)
    626. ) +
    627. scale_fill_manual(values=c("#1B9E77", "#D95F02"))
    628. # display the tip labels of taxa tree
    629. p4 <- p3 + geom_tiplab(size=6, offset=4.0)
    630. # display the LDA of significant OTU.
    631. #p5 <- p4 +
    632. # ggnewscale::new_scale("fill") +
    633. # geom_fruit(
    634. # geom = geom_col,
    635. # mapping = aes(
    636. # x = LDAmean,
    637. # fill = Sign_Group,
    638. # subset = !is.na(LDAmean)
    639. # ),
    640. # orientation = "y",
    641. # offset = 0.3,
    642. # pwidth = 0.5,
    643. # axis.params = list(axis = "x",
    644. # title = "Log10(LDA)",
    645. # title.height = 0.01,
    646. # title.size = 2,
    647. # text.size = 1.8,
    648. # vjust = 1),
    649. # grid.params = list(linetype = 2)
    650. # )
    651. # display the significant (FDR-->pvalue) taxonomy after kruskal.test (default)
    652. #shape = 21,
    653. #scale_size_continuous(range=c(1, 3)) +
    654. p6 <- p4 +
    655. ggnewscale::new_scale("size") +
    656. geom_point(
    657. data=td_filter(!is.na(Sign_Group)),
    658. mapping = aes(size = -log10(pvalue),
    659. fill = Sign_Group,
    660. ),
    661. shape = 21,
    662. ) +
    663. scale_size_continuous(range=c(1, 4)) +
    664. scale_fill_manual(values=c("#1B9E77", "#D95F02"))
    665. svg("diff_analysis_Group2_vs_Group6.svg",width=22, height=22)
    666. #png("differently_expressed_otu.png", width=2000, height=2000)
    667. p6 + theme(
    668. legend.key.height = unit(1.0, "cm"),
    669. legend.key.width = unit(1.0, "cm"),
    670. legend.spacing.y = unit(0.01, "cm"),
    671. legend.text = element_text(size = 20),
    672. legend.title = element_text(size = 20)
    673. #legend.position = c(0.99, 0.01)
    674. )
    675. dev.off()
    676. #-- save the data in an Excel-file --
    677. #install.packages("tidyverse")
    678. #library(tidyverse)
    679. install.packages("dplyr")
    680. install.packages("readr")
    681. library(dplyr)
    682. library(readr)
    683. # Extracting different parts of mpse_abund, assuming these functions are applicable:
    684. abundance_data <- as.data.frame(mp_extract_abundance(mpse_abund))
    685. taxonomy_data <- as.data.frame(mp_extract_taxonomy(mpse_abund))
    686. # Joining the data together by a common identifier, assuming 'OTU' is present in both
    687. combined_data <- left_join(abundance_data, taxonomy_data, by = "OTU")
    688. # Install and load openxlsx package if not already installed
    689. if (!require(openxlsx)) {
    690. install.packages("openxlsx")
    691. library(openxlsx)
    692. }
    693. # Assuming `abundance_data` has lists that need flattening; example for one column:
    694. abundance_data <- abundance_data %>%
    695. mutate(RareAbundanceBySample = sapply(RareAbundanceBySample, function(x) paste(x, collapse = ", "))) %>%
    696. mutate(RareAbundanceByGroup = sapply(RareAbundanceByGroup, function(x) paste(x, collapse = ", ")))
    697. # Check the structure now
    698. str(abundance_data)
    699. # Rename the 'label' column in abundance_data to 'OTU'
    700. colnames(abundance_data)[colnames(abundance_data) == "label"] <- "OTU"
    701. # Merge the cleaned abundance data with taxonomy data
    702. combined_data <- left_join(abundance_data, taxonomy_data, by = "OTU")
    703. # Check the resulting data frame
    704. str(combined_data)
    705. library(openxlsx)
    706. # Write the data frame to an Excel file
    707. write.xlsx(combined_data, "diff_analysis_Group2_vs_Group6.xlsx")
    708. #delete RareAbundanceBySample!
    709. # ------ Deprecated ------
    710. #.data, layout, tree.type, .taxa.class, tiplab.size, offset.abun, pwidth.abun, offset.effsize, pwidth.effsize, group.abun, tiplab.linetype
    711. p <- mpse_abund %>%
    712. mp_plot_diff_res(
    713. group.abun = FALSE,
    714. pwidth.abun=0.01,
    715. offset.abun=0.2,
    716. pwidth.effsize=0.01,
    717. offset.effsize=0.2,
    718. tiplab.size = 5.0,
    719. tiplab.linetype = 3
    720. ) +
    721. scale_fill_manual(values=c("deepskyblue", "orange")) +
    722. scale_fill_manual(
    723. aesthetics = "fill_new", # The fill aes was renamed to "fill_new" for the abundance dotplot layer
    724. values = c("deepskyblue", "orange")
    725. ) +
    726. scale_fill_manual(
    727. aesthetics = "fill_new_new", # The fill aes for hight light layer of tree was renamed to 'fill_new_new'
    728. values = c("#E41A1C", "#377EB8", "#4DAF4A",
    729. "#984EA3", "#FF7F00", "#FFFF33",
    730. "#A65628", "#F781BF", "#00FFFF", "#999999"
    731. )
    732. ) +
    733. theme(
    734. axis.title = element_text(size = 28), # Font size for axis titles
    735. axis.text = element_text(size = 28), # Font size for axis text
    736. plot.title = element_text(size = 16), # Font size for plot title
    737. legend.title = element_text(size = 16), # Font size for legend title
    738. legend.text = element_text(size = 14) # Font size for legend text
    739. )
    740. #p$layers[[2]]$geom <- geom_tiplab(fontsize = 22) # Change 12 to the desired font size
    741. p <- mpse_abund %>%
    742. mp_plot_diff_res(
    743. group.abun = TRUE,
    744. pwidth.abun=0.05,
    745. offset.abun=0.04,
    746. pwidth.effsize=0.3,
    747. offset.effsize=0.6,
    748. tiplab.size = 4.9,
    749. barplot.x = -log10(pvalue)
    750. ) +
    751. scale_fill_manual(values=c("deepskyblue", "orange")) +
    752. scale_fill_manual(
    753. aesthetics = "fill_new", # The fill aes was renamed to "fill_new" for the abundance dotplot layer
    754. values = c("deepskyblue", "orange")
    755. ) +
    756. scale_fill_manual(
    757. aesthetics = "fill_new_new", # The fill aes for hight light layer of tree was renamed to 'fill_new_new'
    758. values = c("#E41A1C", "#377EB8", "#4DAF4A",
    759. "#984EA3", "#FF7F00", "#FFFF33",
    760. "#A65628", "#F781BF", "#00FFFF", "#999999"
    761. )
    762. ) +
    763. theme(
    764. axis.title = element_text(size = 28), # Font size for axis titles
    765. axis.text = element_text(size = 28), # Font size for axis text
    766. plot.title = element_text(size = 28), # Font size for plot title
    767. legend.title = element_text(size = 16), # Font size for legend title
    768. legend.text = element_text(size = 14) # Font size for legend text
    769. )
    770. #png("differently_expressed_otu.png", width=2000, height=2000)
    771. svg("p7.svg",width=22, height=22)
    772. p
    773. dev.off()
    774. p <- mpse_abund %>%
    775. mp_plot_diff_res(
    776. group.abun = FALSE,
    777. pwidth.abun=0.1,
    778. tiplab.size = 0.1,
    779. offset.effsize=0.1,
    780. pwidth.effsize=0.2
    781. ) +
    782. scale_fill_manual(values=c("deepskyblue", "orange")) +
    783. scale_fill_manual(
    784. aesthetics = "fill_new", # The fill aes was renamed to "fill_new" for the abundance dotplot layer
    785. values = c("deepskyblue", "orange")
    786. ) +
    787. scale_fill_manual(
    788. aesthetics = "fill_new_new", # The fill aes for hight light layer of tree was renamed to 'fill_new_new'
    789. values = c("#E41A1C", "#377EB8", "#4DAF4A",
    790. "#984EA3", "#FF7F00", "#FFFF33",
    791. "#A65628", "#F781BF", "#00FFFF", "#999999"
    792. )
    793. )
    794. #svg("p8.svg",width=22, height=22)
    795. png("differently_expressed_otu.png", width=2000, height=2000)
    796. p
    797. dev.off()
    798. f <- mpse_abund %>%
    799. mp_plot_diff_cladogram(
    800. label.size = 2.5,
    801. hilight.alpha = .3,
    802. bg.tree.size = .5,
    803. bg.point.size = 2,
    804. bg.point.stroke = .25
    805. ) +
    806. scale_fill_diff_cladogram( # set the color of different group.
    807. values = c('deepskyblue', 'orange')
    808. ) +
    809. scale_size_continuous(range = c(1, 4))
    810. #png("cladogram.png", width=1000, height=1000)
    811. svg("cladogram.svg", width=10, height=10)
    812. f
    813. dev.off()
    814. ## Extract the OTU table and taxonomy table from the phyloseq object
    815. otu_table <- phyloseq::otu_table(ps.ng.tax_abund) %>% as.data.frame() %>% as.matrix()
    816. tax_table <- phyloseq::tax_table(ps.ng.tax_abund) %>% as.data.frame() %>% as.matrix()
    817. write.csv(otu_table, file="otu_table.csv")
    818. write.csv(tax_table, file="tax_table.csv")
    819. #~/Tools/csv2xls-0.4/csv_to_xls.py otu_table.csv tax_table.csv -d',' -o otu_tax.xls
  2. code of Phyloseq.Rmd

    1. author: ""
    2. date: '`r format(Sys.time(), "%d %m %Y")`'
    3. header-includes:
    4. - \usepackage{color, fancyvrb}
    5. output:
    6. rmdformats::readthedown:
    7. highlight: kate
    8. number_sections : yes
    9. pdf_document:
    10. toc: yes
    11. toc_depth: 2
    12. number_sections : yes
    13. ---
    14. ```{r, echo=FALSE, warning=FALSE}
    15. ## Global options
    16. # TODO: reproduce the html with the additional figure/SVN-files for editing.
    17. # IMPORTANT NOTE: needs before "mkdir figures"
    18. #rmarkdown::render('Phyloseq.Rmd',output_file='Phyloseq.html')
    19. ```
    20. ```{r load-packages, include=FALSE}
    21. library(knitr)
    22. library(rmdformats)
    23. library(readxl)
    24. library(dplyr)
    25. library(kableExtra)
    26. options(max.print="75")
    27. knitr::opts_chunk$set(fig.width=8,
    28. fig.height=6,
    29. eval=TRUE,
    30. cache=TRUE,
    31. echo=TRUE,
    32. prompt=FALSE,
    33. tidy=TRUE,
    34. comment=NA,
    35. message=FALSE,
    36. warning=FALSE)
    37. opts_knit$set(width=85)
    38. # Phyloseq R library
    39. #* Phyloseq web site : https://joey711.github.io/phyloseq/index.html
    40. #* See in particular tutorials for
    41. # - importing data: https://joey711.github.io/phyloseq/import-data.html
    42. # - heat maps: https://joey711.github.io/phyloseq/plot_heatmap-examples.html
    43. ```
    44. # Data
    45. Import raw data and assign sample key:
    46. ```{r, echo=TRUE, warning=FALSE}
    47. #extend map_corrected.txt with Diet and Flora
    48. #setwd("~/DATA/Data_Laura_16S_2/core_diversity_e4753")
    49. map_corrected <- read.csv("map_corrected.txt", sep="\t", row.names=1)
    50. knitr::kable(map_corrected) %>% kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"))
    51. ```
    52. # Prerequisites to be installed
    53. * R : https://pbil.univ-lyon1.fr/CRAN/
    54. * R studio : https://www.rstudio.com/products/rstudio/download/#download
    55. ```R
    56. install.packages("dplyr") # To manipulate dataframes
    57. install.packages("readxl") # To read Excel files into R
    58. install.packages("ggplot2") # for high quality graphics
    59. install.packages("heatmaply")
    60. source("https://bioconductor.org/biocLite.R")
    61. biocLite("phyloseq")
    62. ```
    63. ```{r libraries, echo=TRUE, message=FALSE}
    64. library("readxl") # necessary to import the data from Excel file
    65. library("ggplot2") # graphics
    66. library("picante")
    67. library("microbiome") # data analysis and visualisation
    68. library("phyloseq") # also the basis of data object. Data analysis and visualisation
    69. library("ggpubr") # publication quality figures, based on ggplot2
    70. library("dplyr") # data handling, filter and reformat data frames
    71. library("RColorBrewer") # nice color options
    72. library("heatmaply")
    73. library(vegan)
    74. library(gplots)
    75. ```
    76. # Read the data and create phyloseq objects
    77. Three tables are needed
    78. * OTU
    79. * Taxonomy
    80. * Samples
    81. ```{r, echo=TRUE, warning=FALSE}
    82. #Change your working directory to where the files are located
    83. ps.ng.tax <- import_biom("./table_even42369.biom", "./clustering/rep_set.tre")
    84. sample <- read.csv("./map_corrected.txt", sep="\t", row.names=1)
    85. SAM = sample_data(sample, errorIfNULL = T)
    86. rownames(SAM) <-
    87. c("1","2","3","5","6","7","8","9","10","12","13","14","15","16","17","18","19","20","21","22","23","24","25","26","27","28","29","30","31","32","33","34","35","36","37","38","39","40","41","42","43","44","46","47","48","49","50","51","52","53","55")
    88. ps.ng.tax <- merge_phyloseq(ps.ng.tax, SAM)
    89. print(ps.ng.tax)
    90. colnames(tax_table(ps.ng.tax)) <- c("Domain","Phylum","Class","Order","Family","Genus","Species")
    91. saveRDS(ps.ng.tax, "./ps.ng.tax.rds")
    92. ```
    93. Visualize data
    94. ```{r, echo=TRUE, warning=FALSE}
    95. sample_names(ps.ng.tax)
    96. rank_names(ps.ng.tax)
    97. sample_variables(ps.ng.tax)
    98. ```
    99. Normalize number of reads in each sample using median sequencing depth.
    100. ```{r, echo=TRUE, warning=FALSE}
    101. # RAREFACTION
    102. #set.seed(9242) # This will help in reproducing the filtering and nomalisation.
    103. #ps.ng.tax <- rarefy_even_depth(ps.ng.tax, sample.size = 42369)
    104. #total <- 42369
    105. # NORMALIZE number of reads in each sample using median sequencing depth.
    106. total = median(sample_sums(ps.ng.tax))
    107. #> total
    108. #[1] 42369
    109. standf = function(x, t=total) round(t * (x / sum(x)))
    110. ps.ng.tax = transform_sample_counts(ps.ng.tax, standf)
    111. ps.ng.tax_rel <- microbiome::transform(ps.ng.tax, "compositional")
    112. saveRDS(ps.ng.tax, "./ps.ng.tax.rds")
    113. hmp.meta <- meta(ps.ng.tax)
    114. hmp.meta$sam_name <- rownames(hmp.meta)
    115. ```
    116. # Heatmaps
    117. ```{r, echo=TRUE, warning=FALSE}
    118. #MOVE_FROM_ABOVE: The number of reads used for normalization is **`r sprintf("%.0f", total)`**.
    119. #A basic heatmap using the default parameters.
    120. # plot_heatmap(ps.ng.tax, method = "NMDS", distance = "bray")
    121. #NOTE that giving the correct OTU numbers in the text (1%, 0.5%, ...)!!!
    122. ```
    123. We consider the most abundant OTUs for heatmaps. For example one can only take OTUs that represent at least 1% of reads in at least one sample. Remember we normalized all the sampples to median number of reads (total). We are left with only 168 OTUS which makes the reading much more easy.
    124. ```{r, echo=TRUE, warning=FALSE}
    125. # Custom function to plot a heatmap with the specified sample order
    126. #plot_heatmap_custom <- function(ps, sample_order, method = "NMDS", distance = "bray") {
    127. ps.ng.tax_abund <- phyloseq::filter_taxa(ps.ng.tax, function(x) sum(x > total*0.01) > 0, TRUE)
    128. kable(otu_table(ps.ng.tax_abund)) %>%
    129. kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"))
    130. # Calculate the relative abundance for each sample
    131. ps.ng.tax_abund_rel <- transform_sample_counts(ps.ng.tax_abund, function(x) x / sum(x))
    132. datamat_ = as.data.frame(otu_table(ps.ng.tax_abund))
    133. #datamat <- datamat_[c("1","2","5","6","7", "8","9","10","12","13","14", "15","16","17","18","19","20", "21","22","23","24","25","26","27","28", "29","30","31","32", "33","34","35","36","37","38","39","51", "40","41","42","43","44","46", "47","48","49","50","52","53","55")]
    134. datamat <- datamat_[c("8","9","10","12","13","14", "21","22","23","24","25","26","27","28", "33","34","35","36","37","38","39","51", "47","48","49","50","52","53","55")]
    135. hr <- hclust(as.dist(1-cor(t(datamat), method="pearson")), method="complete")
    136. hc <- hclust(as.dist(1-cor(datamat, method="spearman")), method="complete")
    137. mycl = cutree(hr, h=max(hr$height)/1.08)
    138. mycol = c("YELLOW", "DARKBLUE", "DARKORANGE", "DARKMAGENTA", "DARKCYAN", "DARKRED", "MAROON", "DARKGREEN", "LIGHTBLUE", "PINK", "MAGENTA", "LIGHTCYAN","LIGHTGREEN", "BLUE", "ORANGE", "CYAN", "RED", "GREEN");
    139. mycol = mycol[as.vector(mycl)]
    140. sampleCols <- rep('GREY',ncol(datamat))
    141. #names(sampleCols) <- c("Group1", "Group1", "Group1", "Group1", "Group1", "Group2", "Group2", "Group2", "Group2", "Group2","Group2", "Group3", "Group3", "Group3", "Group3", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1", "Group1")
    142. #sampleCols[colnames(datamat)=='1'] <- '#a6cee3'
    143. #sampleCols[colnames(datamat)=='2'] <- '#a6cee3'
    144. #sampleCols[colnames(datamat)=='5'] <- '#a6cee3'
    145. #sampleCols[colnames(datamat)=='6'] <- '#a6cee3'
    146. #sampleCols[colnames(datamat)=='7'] <- '#a6cee3'
    147. sampleCols[colnames(datamat)=='8'] <- '#1f78b4'
    148. sampleCols[colnames(datamat)=='9'] <- '#1f78b4'
    149. sampleCols[colnames(datamat)=='10'] <- '#1f78b4'
    150. sampleCols[colnames(datamat)=='12'] <- '#1f78b4'
    151. sampleCols[colnames(datamat)=='13'] <- '#1f78b4'
    152. sampleCols[colnames(datamat)=='14'] <- '#1f78b4'
    153. #sampleCols[colnames(datamat)=='15'] <- '#b2df8a'
    154. #sampleCols[colnames(datamat)=='16'] <- '#b2df8a'
    155. #sampleCols[colnames(datamat)=='17'] <- '#b2df8a'
    156. #sampleCols[colnames(datamat)=='18'] <- '#b2df8a'
    157. #sampleCols[colnames(datamat)=='19'] <- '#b2df8a'
    158. #sampleCols[colnames(datamat)=='20'] <- '#b2df8a'
    159. sampleCols[colnames(datamat)=='21'] <- '#33a02c'
    160. sampleCols[colnames(datamat)=='22'] <- '#33a02c'
    161. sampleCols[colnames(datamat)=='23'] <- '#33a02c'
    162. sampleCols[colnames(datamat)=='24'] <- '#33a02c'
    163. sampleCols[colnames(datamat)=='25'] <- '#33a02c'
    164. sampleCols[colnames(datamat)=='26'] <- '#33a02c'
    165. sampleCols[colnames(datamat)=='27'] <- '#33a02c'
    166. sampleCols[colnames(datamat)=='28'] <- '#33a02c'
    167. #sampleCols[colnames(datamat)=='29'] <- '#fb9a99'
    168. #sampleCols[colnames(datamat)=='30'] <- '#fb9a99'
    169. #sampleCols[colnames(datamat)=='31'] <- '#fb9a99'
    170. #sampleCols[colnames(datamat)=='32'] <- '#fb9a99'
    171. sampleCols[colnames(datamat)=='33'] <- '#e31a1c'
    172. sampleCols[colnames(datamat)=='34'] <- '#e31a1c'
    173. sampleCols[colnames(datamat)=='35'] <- '#e31a1c'
    174. sampleCols[colnames(datamat)=='36'] <- '#e31a1c'
    175. sampleCols[colnames(datamat)=='37'] <- '#e31a1c'
    176. sampleCols[colnames(datamat)=='38'] <- '#e31a1c'
    177. sampleCols[colnames(datamat)=='39'] <- '#e31a1c'
    178. sampleCols[colnames(datamat)=='51'] <- '#e31a1c'
    179. #sampleCols[colnames(datamat)=='40'] <- '#cab2d6'
    180. #sampleCols[colnames(datamat)=='41'] <- '#cab2d6'
    181. #sampleCols[colnames(datamat)=='42'] <- '#cab2d6'
    182. #sampleCols[colnames(datamat)=='43'] <- '#cab2d6'
    183. #sampleCols[colnames(datamat)=='44'] <- '#cab2d6'
    184. #sampleCols[colnames(datamat)=='46'] <- '#cab2d6'
    185. sampleCols[colnames(datamat)=='47'] <- '#6a3d9a'
    186. sampleCols[colnames(datamat)=='48'] <- '#6a3d9a'
    187. sampleCols[colnames(datamat)=='49'] <- '#6a3d9a'
    188. sampleCols[colnames(datamat)=='50'] <- '#6a3d9a'
    189. sampleCols[colnames(datamat)=='52'] <- '#6a3d9a'
    190. sampleCols[colnames(datamat)=='53'] <- '#6a3d9a'
    191. sampleCols[colnames(datamat)=='55'] <- '#6a3d9a'
    192. #bluered(75)
    193. #color_pattern <- colorRampPalette(c("blue", "white", "red"))(100)
    194. library(RColorBrewer)
    195. custom_palette <- colorRampPalette(brewer.pal(9, "Blues"))
    196. heatmap_colors <- custom_palette(100)
    197. #colors <- heatmap_color_default(100)
    198. png("figures/heatmap.png", width=1200, height=2400)
    199. #par(mar=c(2, 2, 2, 2)) , lwid=1 lhei=c(0.7, 10)) # Adjust height of color keys keysize=0.3,
    200. heatmap.2(as.matrix(datamat),Rowv=as.dendrogram(hr),Colv = NA, dendrogram = 'row',
    201. scale='row',trace='none',col=heatmap_colors, cexRow=1.2, cexCol=1.5,
    202. RowSideColors = mycol, ColSideColors = sampleCols, srtCol=15, labRow=row.names(datamat), key=TRUE, margins=c(10, 15), lhei=c(0.7, 15), lwid=c(1,8))
    203. dev.off()
    204. ```
    205. ```{r, echo=TRUE, warning=FALSE, fig.cap="Heatmap", out.width = '100%', fig.align= "center"}
    206. knitr::include_graphics("./figures/heatmap.png")
    207. ```
    208. ```{r, echo=FALSE, warning=FALSE}
    209. #It is possible to use different distances and different multivaraite methods. Many different built-in distances can be used.
    210. #dist_methods <- unlist(distanceMethodList)
    211. #print(dist_methods)
    212. ```
    213. \pagebreak
    214. # Taxonomic summary
    215. ## Bar plots in phylum level
    216. ```{r, echo=FALSE, warning=FALSE}
    217. #Make the bargraph nicer by removing OTUs boundaries. This is done by adding ggplot2 modifier.
    218. # 1: uniform color. Color is for the border, fill is for the inside
    219. #ggplot(mtcars, aes(x=as.factor(cyl) )) +
    220. # geom_bar(color="blue", fill=rgb(0.1,0.4,0.5,0.7) )
    221. # 2: Using Hue
    222. #ggplot(mtcars, aes(x=as.factor(cyl), fill=as.factor(cyl) )) +
    223. # geom_bar( ) +
    224. # scale_fill_hue(c = 40) +
    225. # theme(legend.position="none")
    226. # 3: Using RColorBrewer
    227. #ggplot(mtcars, aes(x=as.factor(cyl), fill=as.factor(cyl) )) +
    228. # geom_bar( ) +
    229. # scale_fill_brewer(palette = "Set1") +
    230. # theme(legend.position="none")
    231. # 4: Using greyscale:
    232. #ggplot(mtcars, aes(x=as.factor(cyl), fill=as.factor(cyl) )) +
    233. # geom_bar( ) +
    234. # scale_fill_grey(start = 0.25, end = 0.75) +
    235. # theme(legend.position="none")
    236. # 5: Set manualy
    237. #ggplot(mtcars, aes(x=as.factor(cyl), fill=as.factor(cyl) )) +
    238. # geom_bar( ) +
    239. # scale_fill_manual(values = c("red", "green", "blue") ) +
    240. # theme(legend.position="none")
    241. #NOT SUCCESSFUL!
    242. #allGroupsColors<- c(
    243. # "grey0", "grey50", "dodgerblu", "deepskyblue",
    244. # "red", "darkred", "green", "green4")
    245. # plot_bar(ps.ng.tax_rel, fill="Phylum") +
    246. # geom_bar(stat="identity", position="stack") + scale_color_manual(values = allGroupsColors) #, fill=Phylum + scale_fill_brewer(palette = "Set1")
    247. # ##### Keep only the most abundant phyla and
    248. #ps.ng.tax_most <- subset_taxa(ps.ng.tax_rel, Phylum %in% c("D_1__Actinobacteria")) #1.57
    249. #ps.ng.tax_most <- subset_taxa(ps.ng.tax_rel, Phylum %in% c("D_1__Bacteroidetes")) #27.27436
    250. #ps.ng.tax_most <- subset_taxa(ps.ng.tax_rel, Phylum %in% c("D_1__Cyanobacteria")) #0.02244249
    251. #ps.ng.tax_most <- subset_taxa(ps.ng.tax_rel, Phylum %in% c("D_1__Epsilonbacteraeota")) #0.01309145
    252. #ps.ng.tax_most <- subset_taxa(ps.ng.tax_rel, Phylum %in% c("D_1__Euryarchaeota")) #0.1210024
    253. #ps.ng.tax_most <- subset_taxa(ps.ng.tax_rel, Phylum %in% c("D_1__Firmicutes")) #32.50589
    254. #ps.ng.tax_most <- subset_taxa(ps.ng.tax_rel, Phylum %in% c("D_1__Lentisphaerae")) #0.0001870208
    255. #ps.ng.tax_most <- subset_taxa(ps.ng.tax_rel, Phylum %in% c("D_1__Patescibacteria")) #0.008789976
    256. #ps.ng.tax_most <- subset_taxa(ps.ng.tax_rel, Phylum %in% c("D_1__Planctomycetes")) #0.01365252
    257. #ps.ng.tax_most <- subset_taxa(ps.ng.tax_rel, Phylum %in% c("D_1__Proteobacteria")) #6.769216
    258. #ps.ng.tax_most <- subset_taxa(ps.ng.tax_rel, Phylum %in% c("D_1__Synergistetes")) #0.005049561
    259. #ps.ng.tax_most <- subset_taxa(ps.ng.tax_rel, Phylum %in% c("D_1__Tenericutes")) #0.0005610623
    260. #ps.ng.tax_most <- subset_taxa(ps.ng.tax_rel, Phylum %in% c("D_1__Verrucomicrobia")) #2.076304
    261. #ps.ng.tax_most <- subset_taxa(ps.ng.tax_rel, Phylum %in% c(NA)) #sum(otu_table(ps.ng.tax_most)) = 2.619413
    262. ```
    263. ```{r, echo=TRUE, warning=FALSE}
    264. library(ggplot2)
    265. geom.text.size = 6
    266. theme.size = 8 #(14/5) * geom.text.size
    267. #ps.ng.tax_most <- subset_taxa(ps.ng.tax_rel, Phylum %in% c("D_1__Actinobacteria", "D_1__Bacteroidetes", "D_1__Firmicutes", "D_1__Proteobacteria", "D_1__Verrucomicrobia", NA))
    268. ps.ng.tax_most = phyloseq::filter_taxa(ps.ng.tax_rel, function(x) mean(x) > 0.001, TRUE)
    269. #CONSOLE(OPTIONAL): for sampleid in 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73; do
    270. #echo "otu_table(ps.ng.tax_most)[,${sampleid}]=otu_table(ps.ng.tax_most)[,${sampleid}]/sum(otu_table(ps.ng.tax_most)[,${sampleid}])" done
    271. #OR
    272. ps.ng.tax_most_ = transform_sample_counts(ps.ng.tax_most, function(x) x / sum(x))
    273. ```
    274. ```{r, echo=FALSE, warning=FALSE}
    275. ##--Creating 100% stacked bar plots with less abundant taxa in a sub-category #901--
    276. ##https://github.com/joey711/phyloseq/issues/901
    277. ##ps.ng.tax_most_df <- psmelt(ps.ng.tax_most_) #5986x19
    278. #glom <- tax_glom(ps.ng.tax_most_, taxrank = 'Phylum')
    279. #tax_table(glom) # should list # taxa as # phyla
    280. #data <- psmelt(glom) # create dataframe from phyloseq object
    281. #data$Phylum <- as.character(data$Phylum) #convert to character
    282. ##simple way to rename phyla with < 1% abundance
    283. #data$Phylum[data$Abundance < 0.001] <- "< 0.1% abund."
    284. #
    285. #library(plyr)
    286. #medians <- ddply(data, ~Phylum, function(x) c(median=median(x$Abundance)))
    287. #remainder <- medians[medians$median <= 0.001,]$Phylum
    288. #data[data$Phylum %in% remainder,]$Phylum <- "Phyla < 0.1% abund."
    289. #data$Phylum[data$Abundance < 0.001] <- "Phyla < 0.1% abund."
    290. ##--> data are not used!
    291. #
    292. ##in class level
    293. #glom <- tax_glom(ps.ng.tax_most_, taxrank = 'Class')
    294. #tax_table(glom) # should list # taxa as # phyla
    295. #data <- psmelt(glom) # create dataframe from phyloseq object
    296. #data$Class <- as.character(data$Class) #convert to character
    297. #
    298. ##simple way to rename phyla with < 1% abundance
    299. #data$Class[data$Abundance < 0.001] <- "< 0.1% abund."
    300. #Count = length(unique(data$Class))
    301. #
    302. ##unique(data$Class)
    303. ##data$Class <- factor(data$Class, levels = c("Bacilli", "Bacteroidia", "Verrucomicrobiae", "Clostridia", "Gammaproteobacteria", "Alphaproteobacteria", "Actinobacteria", "Negativicutes", "Erysipelotrichia", "Methanobacteria", "< 0.1% abund."))
    304. ##------- Creating 100% stacked bar plots END --------
    305. library(stringr)
    306. #FITTING1:
    307. # tax_table(ps.ng.tax_most_)[1,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[1,"Domain"], "__")[[1]][2]
    308. # ... ...
    309. # tax_table(ps.ng.tax_most_)[167,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[167,"Species"], "__")[[1]][2]
    310. #ps.ng.tax_most_
    311. #in total [ 89 taxa and 55 samples ]
    312. #otu_table() OTU Table: [ 166 taxa and 54 samples ]
    313. #otu_table() OTU Table: [ 168 taxa and 50 samples ]
    314. #for id in 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166; do
    315. #echo "tax_table(ps.ng.tax_most_)[${id},\"Domain\"] <- str_split(tax_table(ps.ng.tax_most_)[${id},\"Domain\"], \"__\")[[1]][2]"
    316. #echo "tax_table(ps.ng.tax_most_)[${id},\"Phylum\"] <- str_split(tax_table(ps.ng.tax_most_)[${id},\"Phylum\"], \"__\")[[1]][2]"
    317. #echo "tax_table(ps.ng.tax_most_)[${id},\"Class\"] <- str_split(tax_table(ps.ng.tax_most_)[${id},\"Class\"], \"__\")[[1]][2]"
    318. #echo "tax_table(ps.ng.tax_most_)[${id},\"Order\"] <- str_split(tax_table(ps.ng.tax_most_)[${id},\"Order\"], \"__\")[[1]][2]"
    319. #echo "tax_table(ps.ng.tax_most_)[${id},\"Family\"] <- str_split(tax_table(ps.ng.tax_most_)[${id},\"Family\"], \"__\")[[1]][2]"
    320. #echo "tax_table(ps.ng.tax_most_)[${id},\"Genus\"] <- str_split(tax_table(ps.ng.tax_most_)[${id},\"Genus\"], \"__\")[[1]][2]"
    321. #echo "tax_table(ps.ng.tax_most_)[${id},\"Species\"] <- str_split(tax_table(ps.ng.tax_most_)[${id},\"Species\"], \"__\")[[1]][2]"
    322. #done
    323. tax_table(ps.ng.tax_most_)[1,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[1,"Domain"], "__")[[1]][2]
    324. tax_table(ps.ng.tax_most_)[1,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[1,"Phylum"], "__")[[1]][2]
    325. tax_table(ps.ng.tax_most_)[1,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[1,"Class"], "__")[[1]][2]
    326. tax_table(ps.ng.tax_most_)[1,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[1,"Order"], "__")[[1]][2]
    327. tax_table(ps.ng.tax_most_)[1,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[1,"Family"], "__")[[1]][2]
    328. tax_table(ps.ng.tax_most_)[1,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[1,"Genus"], "__")[[1]][2]
    329. tax_table(ps.ng.tax_most_)[1,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[1,"Species"], "__")[[1]][2]
    330. tax_table(ps.ng.tax_most_)[2,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[2,"Domain"], "__")[[1]][2]
    331. tax_table(ps.ng.tax_most_)[2,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[2,"Phylum"], "__")[[1]][2]
    332. tax_table(ps.ng.tax_most_)[2,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[2,"Class"], "__")[[1]][2]
    333. tax_table(ps.ng.tax_most_)[2,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[2,"Order"], "__")[[1]][2]
    334. tax_table(ps.ng.tax_most_)[2,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[2,"Family"], "__")[[1]][2]
    335. tax_table(ps.ng.tax_most_)[2,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[2,"Genus"], "__")[[1]][2]
    336. tax_table(ps.ng.tax_most_)[2,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[2,"Species"], "__")[[1]][2]
    337. tax_table(ps.ng.tax_most_)[3,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[3,"Domain"], "__")[[1]][2]
    338. tax_table(ps.ng.tax_most_)[3,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[3,"Phylum"], "__")[[1]][2]
    339. tax_table(ps.ng.tax_most_)[3,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[3,"Class"], "__")[[1]][2]
    340. tax_table(ps.ng.tax_most_)[3,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[3,"Order"], "__")[[1]][2]
    341. tax_table(ps.ng.tax_most_)[3,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[3,"Family"], "__")[[1]][2]
    342. tax_table(ps.ng.tax_most_)[3,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[3,"Genus"], "__")[[1]][2]
    343. tax_table(ps.ng.tax_most_)[3,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[3,"Species"], "__")[[1]][2]
    344. tax_table(ps.ng.tax_most_)[4,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[4,"Domain"], "__")[[1]][2]
    345. tax_table(ps.ng.tax_most_)[4,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[4,"Phylum"], "__")[[1]][2]
    346. tax_table(ps.ng.tax_most_)[4,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[4,"Class"], "__")[[1]][2]
    347. tax_table(ps.ng.tax_most_)[4,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[4,"Order"], "__")[[1]][2]
    348. tax_table(ps.ng.tax_most_)[4,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[4,"Family"], "__")[[1]][2]
    349. tax_table(ps.ng.tax_most_)[4,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[4,"Genus"], "__")[[1]][2]
    350. tax_table(ps.ng.tax_most_)[4,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[4,"Species"], "__")[[1]][2]
    351. tax_table(ps.ng.tax_most_)[5,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[5,"Domain"], "__")[[1]][2]
    352. tax_table(ps.ng.tax_most_)[5,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[5,"Phylum"], "__")[[1]][2]
    353. tax_table(ps.ng.tax_most_)[5,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[5,"Class"], "__")[[1]][2]
    354. tax_table(ps.ng.tax_most_)[5,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[5,"Order"], "__")[[1]][2]
    355. tax_table(ps.ng.tax_most_)[5,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[5,"Family"], "__")[[1]][2]
    356. tax_table(ps.ng.tax_most_)[5,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[5,"Genus"], "__")[[1]][2]
    357. tax_table(ps.ng.tax_most_)[5,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[5,"Species"], "__")[[1]][2]
    358. tax_table(ps.ng.tax_most_)[6,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[6,"Domain"], "__")[[1]][2]
    359. tax_table(ps.ng.tax_most_)[6,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[6,"Phylum"], "__")[[1]][2]
    360. tax_table(ps.ng.tax_most_)[6,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[6,"Class"], "__")[[1]][2]
    361. tax_table(ps.ng.tax_most_)[6,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[6,"Order"], "__")[[1]][2]
    362. tax_table(ps.ng.tax_most_)[6,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[6,"Family"], "__")[[1]][2]
    363. tax_table(ps.ng.tax_most_)[6,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[6,"Genus"], "__")[[1]][2]
    364. tax_table(ps.ng.tax_most_)[6,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[6,"Species"], "__")[[1]][2]
    365. tax_table(ps.ng.tax_most_)[7,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[7,"Domain"], "__")[[1]][2]
    366. tax_table(ps.ng.tax_most_)[7,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[7,"Phylum"], "__")[[1]][2]
    367. tax_table(ps.ng.tax_most_)[7,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[7,"Class"], "__")[[1]][2]
    368. tax_table(ps.ng.tax_most_)[7,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[7,"Order"], "__")[[1]][2]
    369. tax_table(ps.ng.tax_most_)[7,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[7,"Family"], "__")[[1]][2]
    370. tax_table(ps.ng.tax_most_)[7,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[7,"Genus"], "__")[[1]][2]
    371. tax_table(ps.ng.tax_most_)[7,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[7,"Species"], "__")[[1]][2]
    372. tax_table(ps.ng.tax_most_)[8,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[8,"Domain"], "__")[[1]][2]
    373. tax_table(ps.ng.tax_most_)[8,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[8,"Phylum"], "__")[[1]][2]
    374. tax_table(ps.ng.tax_most_)[8,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[8,"Class"], "__")[[1]][2]
    375. tax_table(ps.ng.tax_most_)[8,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[8,"Order"], "__")[[1]][2]
    376. tax_table(ps.ng.tax_most_)[8,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[8,"Family"], "__")[[1]][2]
    377. tax_table(ps.ng.tax_most_)[8,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[8,"Genus"], "__")[[1]][2]
    378. tax_table(ps.ng.tax_most_)[8,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[8,"Species"], "__")[[1]][2]
    379. tax_table(ps.ng.tax_most_)[9,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[9,"Domain"], "__")[[1]][2]
    380. tax_table(ps.ng.tax_most_)[9,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[9,"Phylum"], "__")[[1]][2]
    381. tax_table(ps.ng.tax_most_)[9,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[9,"Class"], "__")[[1]][2]
    382. tax_table(ps.ng.tax_most_)[9,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[9,"Order"], "__")[[1]][2]
    383. tax_table(ps.ng.tax_most_)[9,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[9,"Family"], "__")[[1]][2]
    384. tax_table(ps.ng.tax_most_)[9,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[9,"Genus"], "__")[[1]][2]
    385. tax_table(ps.ng.tax_most_)[9,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[9,"Species"], "__")[[1]][2]
    386. tax_table(ps.ng.tax_most_)[10,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[10,"Domain"], "__")[[1]][2]
    387. tax_table(ps.ng.tax_most_)[10,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[10,"Phylum"], "__")[[1]][2]
    388. tax_table(ps.ng.tax_most_)[10,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[10,"Class"], "__")[[1]][2]
    389. tax_table(ps.ng.tax_most_)[10,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[10,"Order"], "__")[[1]][2]
    390. tax_table(ps.ng.tax_most_)[10,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[10,"Family"], "__")[[1]][2]
    391. tax_table(ps.ng.tax_most_)[10,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[10,"Genus"], "__")[[1]][2]
    392. tax_table(ps.ng.tax_most_)[10,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[10,"Species"], "__")[[1]][2]
    393. tax_table(ps.ng.tax_most_)[11,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[11,"Domain"], "__")[[1]][2]
    394. tax_table(ps.ng.tax_most_)[11,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[11,"Phylum"], "__")[[1]][2]
    395. tax_table(ps.ng.tax_most_)[11,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[11,"Class"], "__")[[1]][2]
    396. tax_table(ps.ng.tax_most_)[11,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[11,"Order"], "__")[[1]][2]
    397. tax_table(ps.ng.tax_most_)[11,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[11,"Family"], "__")[[1]][2]
    398. tax_table(ps.ng.tax_most_)[11,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[11,"Genus"], "__")[[1]][2]
    399. tax_table(ps.ng.tax_most_)[11,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[11,"Species"], "__")[[1]][2]
    400. tax_table(ps.ng.tax_most_)[12,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[12,"Domain"], "__")[[1]][2]
    401. tax_table(ps.ng.tax_most_)[12,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[12,"Phylum"], "__")[[1]][2]
    402. tax_table(ps.ng.tax_most_)[12,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[12,"Class"], "__")[[1]][2]
    403. tax_table(ps.ng.tax_most_)[12,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[12,"Order"], "__")[[1]][2]
    404. tax_table(ps.ng.tax_most_)[12,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[12,"Family"], "__")[[1]][2]
    405. tax_table(ps.ng.tax_most_)[12,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[12,"Genus"], "__")[[1]][2]
    406. tax_table(ps.ng.tax_most_)[12,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[12,"Species"], "__")[[1]][2]
    407. tax_table(ps.ng.tax_most_)[13,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[13,"Domain"], "__")[[1]][2]
    408. tax_table(ps.ng.tax_most_)[13,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[13,"Phylum"], "__")[[1]][2]
    409. tax_table(ps.ng.tax_most_)[13,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[13,"Class"], "__")[[1]][2]
    410. tax_table(ps.ng.tax_most_)[13,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[13,"Order"], "__")[[1]][2]
    411. tax_table(ps.ng.tax_most_)[13,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[13,"Family"], "__")[[1]][2]
    412. tax_table(ps.ng.tax_most_)[13,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[13,"Genus"], "__")[[1]][2]
    413. tax_table(ps.ng.tax_most_)[13,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[13,"Species"], "__")[[1]][2]
    414. tax_table(ps.ng.tax_most_)[14,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[14,"Domain"], "__")[[1]][2]
    415. tax_table(ps.ng.tax_most_)[14,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[14,"Phylum"], "__")[[1]][2]
    416. tax_table(ps.ng.tax_most_)[14,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[14,"Class"], "__")[[1]][2]
    417. tax_table(ps.ng.tax_most_)[14,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[14,"Order"], "__")[[1]][2]
    418. tax_table(ps.ng.tax_most_)[14,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[14,"Family"], "__")[[1]][2]
    419. tax_table(ps.ng.tax_most_)[14,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[14,"Genus"], "__")[[1]][2]
    420. tax_table(ps.ng.tax_most_)[14,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[14,"Species"], "__")[[1]][2]
    421. tax_table(ps.ng.tax_most_)[15,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[15,"Domain"], "__")[[1]][2]
    422. tax_table(ps.ng.tax_most_)[15,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[15,"Phylum"], "__")[[1]][2]
    423. tax_table(ps.ng.tax_most_)[15,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[15,"Class"], "__")[[1]][2]
    424. tax_table(ps.ng.tax_most_)[15,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[15,"Order"], "__")[[1]][2]
    425. tax_table(ps.ng.tax_most_)[15,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[15,"Family"], "__")[[1]][2]
    426. tax_table(ps.ng.tax_most_)[15,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[15,"Genus"], "__")[[1]][2]
    427. tax_table(ps.ng.tax_most_)[15,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[15,"Species"], "__")[[1]][2]
    428. tax_table(ps.ng.tax_most_)[16,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[16,"Domain"], "__")[[1]][2]
    429. tax_table(ps.ng.tax_most_)[16,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[16,"Phylum"], "__")[[1]][2]
    430. tax_table(ps.ng.tax_most_)[16,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[16,"Class"], "__")[[1]][2]
    431. tax_table(ps.ng.tax_most_)[16,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[16,"Order"], "__")[[1]][2]
    432. tax_table(ps.ng.tax_most_)[16,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[16,"Family"], "__")[[1]][2]
    433. tax_table(ps.ng.tax_most_)[16,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[16,"Genus"], "__")[[1]][2]
    434. tax_table(ps.ng.tax_most_)[16,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[16,"Species"], "__")[[1]][2]
    435. tax_table(ps.ng.tax_most_)[17,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[17,"Domain"], "__")[[1]][2]
    436. tax_table(ps.ng.tax_most_)[17,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[17,"Phylum"], "__")[[1]][2]
    437. tax_table(ps.ng.tax_most_)[17,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[17,"Class"], "__")[[1]][2]
    438. tax_table(ps.ng.tax_most_)[17,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[17,"Order"], "__")[[1]][2]
    439. tax_table(ps.ng.tax_most_)[17,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[17,"Family"], "__")[[1]][2]
    440. tax_table(ps.ng.tax_most_)[17,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[17,"Genus"], "__")[[1]][2]
    441. tax_table(ps.ng.tax_most_)[17,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[17,"Species"], "__")[[1]][2]
    442. tax_table(ps.ng.tax_most_)[18,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[18,"Domain"], "__")[[1]][2]
    443. tax_table(ps.ng.tax_most_)[18,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[18,"Phylum"], "__")[[1]][2]
    444. tax_table(ps.ng.tax_most_)[18,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[18,"Class"], "__")[[1]][2]
    445. tax_table(ps.ng.tax_most_)[18,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[18,"Order"], "__")[[1]][2]
    446. tax_table(ps.ng.tax_most_)[18,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[18,"Family"], "__")[[1]][2]
    447. tax_table(ps.ng.tax_most_)[18,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[18,"Genus"], "__")[[1]][2]
    448. tax_table(ps.ng.tax_most_)[18,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[18,"Species"], "__")[[1]][2]
    449. tax_table(ps.ng.tax_most_)[19,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[19,"Domain"], "__")[[1]][2]
    450. tax_table(ps.ng.tax_most_)[19,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[19,"Phylum"], "__")[[1]][2]
    451. tax_table(ps.ng.tax_most_)[19,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[19,"Class"], "__")[[1]][2]
    452. tax_table(ps.ng.tax_most_)[19,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[19,"Order"], "__")[[1]][2]
    453. tax_table(ps.ng.tax_most_)[19,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[19,"Family"], "__")[[1]][2]
    454. tax_table(ps.ng.tax_most_)[19,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[19,"Genus"], "__")[[1]][2]
    455. tax_table(ps.ng.tax_most_)[19,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[19,"Species"], "__")[[1]][2]
    456. tax_table(ps.ng.tax_most_)[20,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[20,"Domain"], "__")[[1]][2]
    457. tax_table(ps.ng.tax_most_)[20,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[20,"Phylum"], "__")[[1]][2]
    458. tax_table(ps.ng.tax_most_)[20,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[20,"Class"], "__")[[1]][2]
    459. tax_table(ps.ng.tax_most_)[20,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[20,"Order"], "__")[[1]][2]
    460. tax_table(ps.ng.tax_most_)[20,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[20,"Family"], "__")[[1]][2]
    461. tax_table(ps.ng.tax_most_)[20,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[20,"Genus"], "__")[[1]][2]
    462. tax_table(ps.ng.tax_most_)[20,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[20,"Species"], "__")[[1]][2]
    463. tax_table(ps.ng.tax_most_)[21,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[21,"Domain"], "__")[[1]][2]
    464. tax_table(ps.ng.tax_most_)[21,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[21,"Phylum"], "__")[[1]][2]
    465. tax_table(ps.ng.tax_most_)[21,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[21,"Class"], "__")[[1]][2]
    466. tax_table(ps.ng.tax_most_)[21,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[21,"Order"], "__")[[1]][2]
    467. tax_table(ps.ng.tax_most_)[21,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[21,"Family"], "__")[[1]][2]
    468. tax_table(ps.ng.tax_most_)[21,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[21,"Genus"], "__")[[1]][2]
    469. tax_table(ps.ng.tax_most_)[21,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[21,"Species"], "__")[[1]][2]
    470. tax_table(ps.ng.tax_most_)[22,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[22,"Domain"], "__")[[1]][2]
    471. tax_table(ps.ng.tax_most_)[22,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[22,"Phylum"], "__")[[1]][2]
    472. tax_table(ps.ng.tax_most_)[22,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[22,"Class"], "__")[[1]][2]
    473. tax_table(ps.ng.tax_most_)[22,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[22,"Order"], "__")[[1]][2]
    474. tax_table(ps.ng.tax_most_)[22,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[22,"Family"], "__")[[1]][2]
    475. tax_table(ps.ng.tax_most_)[22,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[22,"Genus"], "__")[[1]][2]
    476. tax_table(ps.ng.tax_most_)[22,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[22,"Species"], "__")[[1]][2]
    477. tax_table(ps.ng.tax_most_)[23,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[23,"Domain"], "__")[[1]][2]
    478. tax_table(ps.ng.tax_most_)[23,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[23,"Phylum"], "__")[[1]][2]
    479. tax_table(ps.ng.tax_most_)[23,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[23,"Class"], "__")[[1]][2]
    480. tax_table(ps.ng.tax_most_)[23,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[23,"Order"], "__")[[1]][2]
    481. tax_table(ps.ng.tax_most_)[23,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[23,"Family"], "__")[[1]][2]
    482. tax_table(ps.ng.tax_most_)[23,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[23,"Genus"], "__")[[1]][2]
    483. tax_table(ps.ng.tax_most_)[23,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[23,"Species"], "__")[[1]][2]
    484. tax_table(ps.ng.tax_most_)[24,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[24,"Domain"], "__")[[1]][2]
    485. tax_table(ps.ng.tax_most_)[24,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[24,"Phylum"], "__")[[1]][2]
    486. tax_table(ps.ng.tax_most_)[24,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[24,"Class"], "__")[[1]][2]
    487. tax_table(ps.ng.tax_most_)[24,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[24,"Order"], "__")[[1]][2]
    488. tax_table(ps.ng.tax_most_)[24,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[24,"Family"], "__")[[1]][2]
    489. tax_table(ps.ng.tax_most_)[24,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[24,"Genus"], "__")[[1]][2]
    490. tax_table(ps.ng.tax_most_)[24,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[24,"Species"], "__")[[1]][2]
    491. tax_table(ps.ng.tax_most_)[25,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[25,"Domain"], "__")[[1]][2]
    492. tax_table(ps.ng.tax_most_)[25,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[25,"Phylum"], "__")[[1]][2]
    493. tax_table(ps.ng.tax_most_)[25,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[25,"Class"], "__")[[1]][2]
    494. tax_table(ps.ng.tax_most_)[25,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[25,"Order"], "__")[[1]][2]
    495. tax_table(ps.ng.tax_most_)[25,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[25,"Family"], "__")[[1]][2]
    496. tax_table(ps.ng.tax_most_)[25,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[25,"Genus"], "__")[[1]][2]
    497. tax_table(ps.ng.tax_most_)[25,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[25,"Species"], "__")[[1]][2]
    498. tax_table(ps.ng.tax_most_)[26,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[26,"Domain"], "__")[[1]][2]
    499. tax_table(ps.ng.tax_most_)[26,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[26,"Phylum"], "__")[[1]][2]
    500. tax_table(ps.ng.tax_most_)[26,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[26,"Class"], "__")[[1]][2]
    501. tax_table(ps.ng.tax_most_)[26,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[26,"Order"], "__")[[1]][2]
    502. tax_table(ps.ng.tax_most_)[26,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[26,"Family"], "__")[[1]][2]
    503. tax_table(ps.ng.tax_most_)[26,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[26,"Genus"], "__")[[1]][2]
    504. tax_table(ps.ng.tax_most_)[26,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[26,"Species"], "__")[[1]][2]
    505. tax_table(ps.ng.tax_most_)[27,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[27,"Domain"], "__")[[1]][2]
    506. tax_table(ps.ng.tax_most_)[27,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[27,"Phylum"], "__")[[1]][2]
    507. tax_table(ps.ng.tax_most_)[27,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[27,"Class"], "__")[[1]][2]
    508. tax_table(ps.ng.tax_most_)[27,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[27,"Order"], "__")[[1]][2]
    509. tax_table(ps.ng.tax_most_)[27,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[27,"Family"], "__")[[1]][2]
    510. tax_table(ps.ng.tax_most_)[27,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[27,"Genus"], "__")[[1]][2]
    511. tax_table(ps.ng.tax_most_)[27,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[27,"Species"], "__")[[1]][2]
    512. tax_table(ps.ng.tax_most_)[28,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[28,"Domain"], "__")[[1]][2]
    513. tax_table(ps.ng.tax_most_)[28,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[28,"Phylum"], "__")[[1]][2]
    514. tax_table(ps.ng.tax_most_)[28,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[28,"Class"], "__")[[1]][2]
    515. tax_table(ps.ng.tax_most_)[28,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[28,"Order"], "__")[[1]][2]
    516. tax_table(ps.ng.tax_most_)[28,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[28,"Family"], "__")[[1]][2]
    517. tax_table(ps.ng.tax_most_)[28,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[28,"Genus"], "__")[[1]][2]
    518. tax_table(ps.ng.tax_most_)[28,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[28,"Species"], "__")[[1]][2]
    519. tax_table(ps.ng.tax_most_)[29,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[29,"Domain"], "__")[[1]][2]
    520. tax_table(ps.ng.tax_most_)[29,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[29,"Phylum"], "__")[[1]][2]
    521. tax_table(ps.ng.tax_most_)[29,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[29,"Class"], "__")[[1]][2]
    522. tax_table(ps.ng.tax_most_)[29,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[29,"Order"], "__")[[1]][2]
    523. tax_table(ps.ng.tax_most_)[29,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[29,"Family"], "__")[[1]][2]
    524. tax_table(ps.ng.tax_most_)[29,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[29,"Genus"], "__")[[1]][2]
    525. tax_table(ps.ng.tax_most_)[29,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[29,"Species"], "__")[[1]][2]
    526. tax_table(ps.ng.tax_most_)[30,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[30,"Domain"], "__")[[1]][2]
    527. tax_table(ps.ng.tax_most_)[30,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[30,"Phylum"], "__")[[1]][2]
    528. tax_table(ps.ng.tax_most_)[30,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[30,"Class"], "__")[[1]][2]
    529. tax_table(ps.ng.tax_most_)[30,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[30,"Order"], "__")[[1]][2]
    530. tax_table(ps.ng.tax_most_)[30,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[30,"Family"], "__")[[1]][2]
    531. tax_table(ps.ng.tax_most_)[30,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[30,"Genus"], "__")[[1]][2]
    532. tax_table(ps.ng.tax_most_)[30,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[30,"Species"], "__")[[1]][2]
    533. tax_table(ps.ng.tax_most_)[31,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[31,"Domain"], "__")[[1]][2]
    534. tax_table(ps.ng.tax_most_)[31,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[31,"Phylum"], "__")[[1]][2]
    535. tax_table(ps.ng.tax_most_)[31,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[31,"Class"], "__")[[1]][2]
    536. tax_table(ps.ng.tax_most_)[31,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[31,"Order"], "__")[[1]][2]
    537. tax_table(ps.ng.tax_most_)[31,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[31,"Family"], "__")[[1]][2]
    538. tax_table(ps.ng.tax_most_)[31,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[31,"Genus"], "__")[[1]][2]
    539. tax_table(ps.ng.tax_most_)[31,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[31,"Species"], "__")[[1]][2]
    540. tax_table(ps.ng.tax_most_)[32,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[32,"Domain"], "__")[[1]][2]
    541. tax_table(ps.ng.tax_most_)[32,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[32,"Phylum"], "__")[[1]][2]
    542. tax_table(ps.ng.tax_most_)[32,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[32,"Class"], "__")[[1]][2]
    543. tax_table(ps.ng.tax_most_)[32,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[32,"Order"], "__")[[1]][2]
    544. tax_table(ps.ng.tax_most_)[32,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[32,"Family"], "__")[[1]][2]
    545. tax_table(ps.ng.tax_most_)[32,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[32,"Genus"], "__")[[1]][2]
    546. tax_table(ps.ng.tax_most_)[32,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[32,"Species"], "__")[[1]][2]
    547. tax_table(ps.ng.tax_most_)[33,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[33,"Domain"], "__")[[1]][2]
    548. tax_table(ps.ng.tax_most_)[33,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[33,"Phylum"], "__")[[1]][2]
    549. tax_table(ps.ng.tax_most_)[33,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[33,"Class"], "__")[[1]][2]
    550. tax_table(ps.ng.tax_most_)[33,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[33,"Order"], "__")[[1]][2]
    551. tax_table(ps.ng.tax_most_)[33,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[33,"Family"], "__")[[1]][2]
    552. tax_table(ps.ng.tax_most_)[33,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[33,"Genus"], "__")[[1]][2]
    553. tax_table(ps.ng.tax_most_)[33,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[33,"Species"], "__")[[1]][2]
    554. tax_table(ps.ng.tax_most_)[34,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[34,"Domain"], "__")[[1]][2]
    555. tax_table(ps.ng.tax_most_)[34,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[34,"Phylum"], "__")[[1]][2]
    556. tax_table(ps.ng.tax_most_)[34,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[34,"Class"], "__")[[1]][2]
    557. tax_table(ps.ng.tax_most_)[34,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[34,"Order"], "__")[[1]][2]
    558. tax_table(ps.ng.tax_most_)[34,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[34,"Family"], "__")[[1]][2]
    559. tax_table(ps.ng.tax_most_)[34,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[34,"Genus"], "__")[[1]][2]
    560. tax_table(ps.ng.tax_most_)[34,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[34,"Species"], "__")[[1]][2]
    561. tax_table(ps.ng.tax_most_)[35,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[35,"Domain"], "__")[[1]][2]
    562. tax_table(ps.ng.tax_most_)[35,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[35,"Phylum"], "__")[[1]][2]
    563. tax_table(ps.ng.tax_most_)[35,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[35,"Class"], "__")[[1]][2]
    564. tax_table(ps.ng.tax_most_)[35,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[35,"Order"], "__")[[1]][2]
    565. tax_table(ps.ng.tax_most_)[35,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[35,"Family"], "__")[[1]][2]
    566. tax_table(ps.ng.tax_most_)[35,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[35,"Genus"], "__")[[1]][2]
    567. tax_table(ps.ng.tax_most_)[35,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[35,"Species"], "__")[[1]][2]
    568. tax_table(ps.ng.tax_most_)[36,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[36,"Domain"], "__")[[1]][2]
    569. tax_table(ps.ng.tax_most_)[36,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[36,"Phylum"], "__")[[1]][2]
    570. tax_table(ps.ng.tax_most_)[36,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[36,"Class"], "__")[[1]][2]
    571. tax_table(ps.ng.tax_most_)[36,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[36,"Order"], "__")[[1]][2]
    572. tax_table(ps.ng.tax_most_)[36,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[36,"Family"], "__")[[1]][2]
    573. tax_table(ps.ng.tax_most_)[36,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[36,"Genus"], "__")[[1]][2]
    574. tax_table(ps.ng.tax_most_)[36,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[36,"Species"], "__")[[1]][2]
    575. tax_table(ps.ng.tax_most_)[37,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[37,"Domain"], "__")[[1]][2]
    576. tax_table(ps.ng.tax_most_)[37,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[37,"Phylum"], "__")[[1]][2]
    577. tax_table(ps.ng.tax_most_)[37,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[37,"Class"], "__")[[1]][2]
    578. tax_table(ps.ng.tax_most_)[37,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[37,"Order"], "__")[[1]][2]
    579. tax_table(ps.ng.tax_most_)[37,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[37,"Family"], "__")[[1]][2]
    580. tax_table(ps.ng.tax_most_)[37,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[37,"Genus"], "__")[[1]][2]
    581. tax_table(ps.ng.tax_most_)[37,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[37,"Species"], "__")[[1]][2]
    582. tax_table(ps.ng.tax_most_)[38,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[38,"Domain"], "__")[[1]][2]
    583. tax_table(ps.ng.tax_most_)[38,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[38,"Phylum"], "__")[[1]][2]
    584. tax_table(ps.ng.tax_most_)[38,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[38,"Class"], "__")[[1]][2]
    585. tax_table(ps.ng.tax_most_)[38,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[38,"Order"], "__")[[1]][2]
    586. tax_table(ps.ng.tax_most_)[38,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[38,"Family"], "__")[[1]][2]
    587. tax_table(ps.ng.tax_most_)[38,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[38,"Genus"], "__")[[1]][2]
    588. tax_table(ps.ng.tax_most_)[38,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[38,"Species"], "__")[[1]][2]
    589. tax_table(ps.ng.tax_most_)[39,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[39,"Domain"], "__")[[1]][2]
    590. tax_table(ps.ng.tax_most_)[39,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[39,"Phylum"], "__")[[1]][2]
    591. tax_table(ps.ng.tax_most_)[39,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[39,"Class"], "__")[[1]][2]
    592. tax_table(ps.ng.tax_most_)[39,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[39,"Order"], "__")[[1]][2]
    593. tax_table(ps.ng.tax_most_)[39,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[39,"Family"], "__")[[1]][2]
    594. tax_table(ps.ng.tax_most_)[39,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[39,"Genus"], "__")[[1]][2]
    595. tax_table(ps.ng.tax_most_)[39,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[39,"Species"], "__")[[1]][2]
    596. tax_table(ps.ng.tax_most_)[40,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[40,"Domain"], "__")[[1]][2]
    597. tax_table(ps.ng.tax_most_)[40,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[40,"Phylum"], "__")[[1]][2]
    598. tax_table(ps.ng.tax_most_)[40,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[40,"Class"], "__")[[1]][2]
    599. tax_table(ps.ng.tax_most_)[40,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[40,"Order"], "__")[[1]][2]
    600. tax_table(ps.ng.tax_most_)[40,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[40,"Family"], "__")[[1]][2]
    601. tax_table(ps.ng.tax_most_)[40,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[40,"Genus"], "__")[[1]][2]
    602. tax_table(ps.ng.tax_most_)[40,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[40,"Species"], "__")[[1]][2]
    603. tax_table(ps.ng.tax_most_)[41,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[41,"Domain"], "__")[[1]][2]
    604. tax_table(ps.ng.tax_most_)[41,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[41,"Phylum"], "__")[[1]][2]
    605. tax_table(ps.ng.tax_most_)[41,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[41,"Class"], "__")[[1]][2]
    606. tax_table(ps.ng.tax_most_)[41,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[41,"Order"], "__")[[1]][2]
    607. tax_table(ps.ng.tax_most_)[41,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[41,"Family"], "__")[[1]][2]
    608. tax_table(ps.ng.tax_most_)[41,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[41,"Genus"], "__")[[1]][2]
    609. tax_table(ps.ng.tax_most_)[41,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[41,"Species"], "__")[[1]][2]
    610. tax_table(ps.ng.tax_most_)[42,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[42,"Domain"], "__")[[1]][2]
    611. tax_table(ps.ng.tax_most_)[42,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[42,"Phylum"], "__")[[1]][2]
    612. tax_table(ps.ng.tax_most_)[42,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[42,"Class"], "__")[[1]][2]
    613. tax_table(ps.ng.tax_most_)[42,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[42,"Order"], "__")[[1]][2]
    614. tax_table(ps.ng.tax_most_)[42,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[42,"Family"], "__")[[1]][2]
    615. tax_table(ps.ng.tax_most_)[42,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[42,"Genus"], "__")[[1]][2]
    616. tax_table(ps.ng.tax_most_)[42,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[42,"Species"], "__")[[1]][2]
    617. tax_table(ps.ng.tax_most_)[43,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[43,"Domain"], "__")[[1]][2]
    618. tax_table(ps.ng.tax_most_)[43,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[43,"Phylum"], "__")[[1]][2]
    619. tax_table(ps.ng.tax_most_)[43,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[43,"Class"], "__")[[1]][2]
    620. tax_table(ps.ng.tax_most_)[43,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[43,"Order"], "__")[[1]][2]
    621. tax_table(ps.ng.tax_most_)[43,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[43,"Family"], "__")[[1]][2]
    622. tax_table(ps.ng.tax_most_)[43,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[43,"Genus"], "__")[[1]][2]
    623. tax_table(ps.ng.tax_most_)[43,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[43,"Species"], "__")[[1]][2]
    624. tax_table(ps.ng.tax_most_)[44,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[44,"Domain"], "__")[[1]][2]
    625. tax_table(ps.ng.tax_most_)[44,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[44,"Phylum"], "__")[[1]][2]
    626. tax_table(ps.ng.tax_most_)[44,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[44,"Class"], "__")[[1]][2]
    627. tax_table(ps.ng.tax_most_)[44,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[44,"Order"], "__")[[1]][2]
    628. tax_table(ps.ng.tax_most_)[44,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[44,"Family"], "__")[[1]][2]
    629. tax_table(ps.ng.tax_most_)[44,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[44,"Genus"], "__")[[1]][2]
    630. tax_table(ps.ng.tax_most_)[44,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[44,"Species"], "__")[[1]][2]
    631. tax_table(ps.ng.tax_most_)[45,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[45,"Domain"], "__")[[1]][2]
    632. tax_table(ps.ng.tax_most_)[45,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[45,"Phylum"], "__")[[1]][2]
    633. tax_table(ps.ng.tax_most_)[45,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[45,"Class"], "__")[[1]][2]
    634. tax_table(ps.ng.tax_most_)[45,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[45,"Order"], "__")[[1]][2]
    635. tax_table(ps.ng.tax_most_)[45,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[45,"Family"], "__")[[1]][2]
    636. tax_table(ps.ng.tax_most_)[45,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[45,"Genus"], "__")[[1]][2]
    637. tax_table(ps.ng.tax_most_)[45,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[45,"Species"], "__")[[1]][2]
    638. tax_table(ps.ng.tax_most_)[46,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[46,"Domain"], "__")[[1]][2]
    639. tax_table(ps.ng.tax_most_)[46,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[46,"Phylum"], "__")[[1]][2]
    640. tax_table(ps.ng.tax_most_)[46,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[46,"Class"], "__")[[1]][2]
    641. tax_table(ps.ng.tax_most_)[46,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[46,"Order"], "__")[[1]][2]
    642. tax_table(ps.ng.tax_most_)[46,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[46,"Family"], "__")[[1]][2]
    643. tax_table(ps.ng.tax_most_)[46,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[46,"Genus"], "__")[[1]][2]
    644. tax_table(ps.ng.tax_most_)[46,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[46,"Species"], "__")[[1]][2]
    645. tax_table(ps.ng.tax_most_)[47,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[47,"Domain"], "__")[[1]][2]
    646. tax_table(ps.ng.tax_most_)[47,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[47,"Phylum"], "__")[[1]][2]
    647. tax_table(ps.ng.tax_most_)[47,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[47,"Class"], "__")[[1]][2]
    648. tax_table(ps.ng.tax_most_)[47,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[47,"Order"], "__")[[1]][2]
    649. tax_table(ps.ng.tax_most_)[47,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[47,"Family"], "__")[[1]][2]
    650. tax_table(ps.ng.tax_most_)[47,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[47,"Genus"], "__")[[1]][2]
    651. tax_table(ps.ng.tax_most_)[47,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[47,"Species"], "__")[[1]][2]
    652. tax_table(ps.ng.tax_most_)[48,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[48,"Domain"], "__")[[1]][2]
    653. tax_table(ps.ng.tax_most_)[48,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[48,"Phylum"], "__")[[1]][2]
    654. tax_table(ps.ng.tax_most_)[48,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[48,"Class"], "__")[[1]][2]
    655. tax_table(ps.ng.tax_most_)[48,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[48,"Order"], "__")[[1]][2]
    656. tax_table(ps.ng.tax_most_)[48,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[48,"Family"], "__")[[1]][2]
    657. tax_table(ps.ng.tax_most_)[48,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[48,"Genus"], "__")[[1]][2]
    658. tax_table(ps.ng.tax_most_)[48,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[48,"Species"], "__")[[1]][2]
    659. tax_table(ps.ng.tax_most_)[49,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[49,"Domain"], "__")[[1]][2]
    660. tax_table(ps.ng.tax_most_)[49,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[49,"Phylum"], "__")[[1]][2]
    661. tax_table(ps.ng.tax_most_)[49,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[49,"Class"], "__")[[1]][2]
    662. tax_table(ps.ng.tax_most_)[49,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[49,"Order"], "__")[[1]][2]
    663. tax_table(ps.ng.tax_most_)[49,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[49,"Family"], "__")[[1]][2]
    664. tax_table(ps.ng.tax_most_)[49,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[49,"Genus"], "__")[[1]][2]
    665. tax_table(ps.ng.tax_most_)[49,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[49,"Species"], "__")[[1]][2]
    666. tax_table(ps.ng.tax_most_)[50,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[50,"Domain"], "__")[[1]][2]
    667. tax_table(ps.ng.tax_most_)[50,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[50,"Phylum"], "__")[[1]][2]
    668. tax_table(ps.ng.tax_most_)[50,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[50,"Class"], "__")[[1]][2]
    669. tax_table(ps.ng.tax_most_)[50,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[50,"Order"], "__")[[1]][2]
    670. tax_table(ps.ng.tax_most_)[50,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[50,"Family"], "__")[[1]][2]
    671. tax_table(ps.ng.tax_most_)[50,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[50,"Genus"], "__")[[1]][2]
    672. tax_table(ps.ng.tax_most_)[50,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[50,"Species"], "__")[[1]][2]
    673. tax_table(ps.ng.tax_most_)[51,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[51,"Domain"], "__")[[1]][2]
    674. tax_table(ps.ng.tax_most_)[51,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[51,"Phylum"], "__")[[1]][2]
    675. tax_table(ps.ng.tax_most_)[51,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[51,"Class"], "__")[[1]][2]
    676. tax_table(ps.ng.tax_most_)[51,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[51,"Order"], "__")[[1]][2]
    677. tax_table(ps.ng.tax_most_)[51,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[51,"Family"], "__")[[1]][2]
    678. tax_table(ps.ng.tax_most_)[51,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[51,"Genus"], "__")[[1]][2]
    679. tax_table(ps.ng.tax_most_)[51,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[51,"Species"], "__")[[1]][2]
    680. tax_table(ps.ng.tax_most_)[52,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[52,"Domain"], "__")[[1]][2]
    681. tax_table(ps.ng.tax_most_)[52,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[52,"Phylum"], "__")[[1]][2]
    682. tax_table(ps.ng.tax_most_)[52,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[52,"Class"], "__")[[1]][2]
    683. tax_table(ps.ng.tax_most_)[52,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[52,"Order"], "__")[[1]][2]
    684. tax_table(ps.ng.tax_most_)[52,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[52,"Family"], "__")[[1]][2]
    685. tax_table(ps.ng.tax_most_)[52,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[52,"Genus"], "__")[[1]][2]
    686. tax_table(ps.ng.tax_most_)[52,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[52,"Species"], "__")[[1]][2]
    687. tax_table(ps.ng.tax_most_)[53,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[53,"Domain"], "__")[[1]][2]
    688. tax_table(ps.ng.tax_most_)[53,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[53,"Phylum"], "__")[[1]][2]
    689. tax_table(ps.ng.tax_most_)[53,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[53,"Class"], "__")[[1]][2]
    690. tax_table(ps.ng.tax_most_)[53,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[53,"Order"], "__")[[1]][2]
    691. tax_table(ps.ng.tax_most_)[53,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[53,"Family"], "__")[[1]][2]
    692. tax_table(ps.ng.tax_most_)[53,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[53,"Genus"], "__")[[1]][2]
    693. tax_table(ps.ng.tax_most_)[53,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[53,"Species"], "__")[[1]][2]
    694. tax_table(ps.ng.tax_most_)[54,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[54,"Domain"], "__")[[1]][2]
    695. tax_table(ps.ng.tax_most_)[54,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[54,"Phylum"], "__")[[1]][2]
    696. tax_table(ps.ng.tax_most_)[54,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[54,"Class"], "__")[[1]][2]
    697. tax_table(ps.ng.tax_most_)[54,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[54,"Order"], "__")[[1]][2]
    698. tax_table(ps.ng.tax_most_)[54,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[54,"Family"], "__")[[1]][2]
    699. tax_table(ps.ng.tax_most_)[54,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[54,"Genus"], "__")[[1]][2]
    700. tax_table(ps.ng.tax_most_)[54,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[54,"Species"], "__")[[1]][2]
    701. tax_table(ps.ng.tax_most_)[55,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[55,"Domain"], "__")[[1]][2]
    702. tax_table(ps.ng.tax_most_)[55,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[55,"Phylum"], "__")[[1]][2]
    703. tax_table(ps.ng.tax_most_)[55,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[55,"Class"], "__")[[1]][2]
    704. tax_table(ps.ng.tax_most_)[55,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[55,"Order"], "__")[[1]][2]
    705. tax_table(ps.ng.tax_most_)[55,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[55,"Family"], "__")[[1]][2]
    706. tax_table(ps.ng.tax_most_)[55,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[55,"Genus"], "__")[[1]][2]
    707. tax_table(ps.ng.tax_most_)[55,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[55,"Species"], "__")[[1]][2]
    708. tax_table(ps.ng.tax_most_)[56,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[56,"Domain"], "__")[[1]][2]
    709. tax_table(ps.ng.tax_most_)[56,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[56,"Phylum"], "__")[[1]][2]
    710. tax_table(ps.ng.tax_most_)[56,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[56,"Class"], "__")[[1]][2]
    711. tax_table(ps.ng.tax_most_)[56,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[56,"Order"], "__")[[1]][2]
    712. tax_table(ps.ng.tax_most_)[56,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[56,"Family"], "__")[[1]][2]
    713. tax_table(ps.ng.tax_most_)[56,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[56,"Genus"], "__")[[1]][2]
    714. tax_table(ps.ng.tax_most_)[56,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[56,"Species"], "__")[[1]][2]
    715. tax_table(ps.ng.tax_most_)[57,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[57,"Domain"], "__")[[1]][2]
    716. tax_table(ps.ng.tax_most_)[57,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[57,"Phylum"], "__")[[1]][2]
    717. tax_table(ps.ng.tax_most_)[57,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[57,"Class"], "__")[[1]][2]
    718. tax_table(ps.ng.tax_most_)[57,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[57,"Order"], "__")[[1]][2]
    719. tax_table(ps.ng.tax_most_)[57,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[57,"Family"], "__")[[1]][2]
    720. tax_table(ps.ng.tax_most_)[57,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[57,"Genus"], "__")[[1]][2]
    721. tax_table(ps.ng.tax_most_)[57,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[57,"Species"], "__")[[1]][2]
    722. tax_table(ps.ng.tax_most_)[58,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[58,"Domain"], "__")[[1]][2]
    723. tax_table(ps.ng.tax_most_)[58,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[58,"Phylum"], "__")[[1]][2]
    724. tax_table(ps.ng.tax_most_)[58,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[58,"Class"], "__")[[1]][2]
    725. tax_table(ps.ng.tax_most_)[58,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[58,"Order"], "__")[[1]][2]
    726. tax_table(ps.ng.tax_most_)[58,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[58,"Family"], "__")[[1]][2]
    727. tax_table(ps.ng.tax_most_)[58,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[58,"Genus"], "__")[[1]][2]
    728. tax_table(ps.ng.tax_most_)[58,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[58,"Species"], "__")[[1]][2]
    729. tax_table(ps.ng.tax_most_)[59,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[59,"Domain"], "__")[[1]][2]
    730. tax_table(ps.ng.tax_most_)[59,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[59,"Phylum"], "__")[[1]][2]
    731. tax_table(ps.ng.tax_most_)[59,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[59,"Class"], "__")[[1]][2]
    732. tax_table(ps.ng.tax_most_)[59,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[59,"Order"], "__")[[1]][2]
    733. tax_table(ps.ng.tax_most_)[59,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[59,"Family"], "__")[[1]][2]
    734. tax_table(ps.ng.tax_most_)[59,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[59,"Genus"], "__")[[1]][2]
    735. tax_table(ps.ng.tax_most_)[59,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[59,"Species"], "__")[[1]][2]
    736. tax_table(ps.ng.tax_most_)[60,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[60,"Domain"], "__")[[1]][2]
    737. tax_table(ps.ng.tax_most_)[60,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[60,"Phylum"], "__")[[1]][2]
    738. tax_table(ps.ng.tax_most_)[60,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[60,"Class"], "__")[[1]][2]
    739. tax_table(ps.ng.tax_most_)[60,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[60,"Order"], "__")[[1]][2]
    740. tax_table(ps.ng.tax_most_)[60,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[60,"Family"], "__")[[1]][2]
    741. tax_table(ps.ng.tax_most_)[60,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[60,"Genus"], "__")[[1]][2]
    742. tax_table(ps.ng.tax_most_)[60,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[60,"Species"], "__")[[1]][2]
    743. tax_table(ps.ng.tax_most_)[61,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[61,"Domain"], "__")[[1]][2]
    744. tax_table(ps.ng.tax_most_)[61,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[61,"Phylum"], "__")[[1]][2]
    745. tax_table(ps.ng.tax_most_)[61,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[61,"Class"], "__")[[1]][2]
    746. tax_table(ps.ng.tax_most_)[61,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[61,"Order"], "__")[[1]][2]
    747. tax_table(ps.ng.tax_most_)[61,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[61,"Family"], "__")[[1]][2]
    748. tax_table(ps.ng.tax_most_)[61,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[61,"Genus"], "__")[[1]][2]
    749. tax_table(ps.ng.tax_most_)[61,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[61,"Species"], "__")[[1]][2]
    750. tax_table(ps.ng.tax_most_)[62,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[62,"Domain"], "__")[[1]][2]
    751. tax_table(ps.ng.tax_most_)[62,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[62,"Phylum"], "__")[[1]][2]
    752. tax_table(ps.ng.tax_most_)[62,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[62,"Class"], "__")[[1]][2]
    753. tax_table(ps.ng.tax_most_)[62,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[62,"Order"], "__")[[1]][2]
    754. tax_table(ps.ng.tax_most_)[62,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[62,"Family"], "__")[[1]][2]
    755. tax_table(ps.ng.tax_most_)[62,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[62,"Genus"], "__")[[1]][2]
    756. tax_table(ps.ng.tax_most_)[62,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[62,"Species"], "__")[[1]][2]
    757. tax_table(ps.ng.tax_most_)[63,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[63,"Domain"], "__")[[1]][2]
    758. tax_table(ps.ng.tax_most_)[63,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[63,"Phylum"], "__")[[1]][2]
    759. tax_table(ps.ng.tax_most_)[63,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[63,"Class"], "__")[[1]][2]
    760. tax_table(ps.ng.tax_most_)[63,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[63,"Order"], "__")[[1]][2]
    761. tax_table(ps.ng.tax_most_)[63,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[63,"Family"], "__")[[1]][2]
    762. tax_table(ps.ng.tax_most_)[63,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[63,"Genus"], "__")[[1]][2]
    763. tax_table(ps.ng.tax_most_)[63,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[63,"Species"], "__")[[1]][2]
    764. tax_table(ps.ng.tax_most_)[64,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[64,"Domain"], "__")[[1]][2]
    765. tax_table(ps.ng.tax_most_)[64,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[64,"Phylum"], "__")[[1]][2]
    766. tax_table(ps.ng.tax_most_)[64,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[64,"Class"], "__")[[1]][2]
    767. tax_table(ps.ng.tax_most_)[64,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[64,"Order"], "__")[[1]][2]
    768. tax_table(ps.ng.tax_most_)[64,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[64,"Family"], "__")[[1]][2]
    769. tax_table(ps.ng.tax_most_)[64,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[64,"Genus"], "__")[[1]][2]
    770. tax_table(ps.ng.tax_most_)[64,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[64,"Species"], "__")[[1]][2]
    771. tax_table(ps.ng.tax_most_)[65,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[65,"Domain"], "__")[[1]][2]
    772. tax_table(ps.ng.tax_most_)[65,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[65,"Phylum"], "__")[[1]][2]
    773. tax_table(ps.ng.tax_most_)[65,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[65,"Class"], "__")[[1]][2]
    774. tax_table(ps.ng.tax_most_)[65,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[65,"Order"], "__")[[1]][2]
    775. tax_table(ps.ng.tax_most_)[65,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[65,"Family"], "__")[[1]][2]
    776. tax_table(ps.ng.tax_most_)[65,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[65,"Genus"], "__")[[1]][2]
    777. tax_table(ps.ng.tax_most_)[65,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[65,"Species"], "__")[[1]][2]
    778. tax_table(ps.ng.tax_most_)[66,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[66,"Domain"], "__")[[1]][2]
    779. tax_table(ps.ng.tax_most_)[66,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[66,"Phylum"], "__")[[1]][2]
    780. tax_table(ps.ng.tax_most_)[66,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[66,"Class"], "__")[[1]][2]
    781. tax_table(ps.ng.tax_most_)[66,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[66,"Order"], "__")[[1]][2]
    782. tax_table(ps.ng.tax_most_)[66,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[66,"Family"], "__")[[1]][2]
    783. tax_table(ps.ng.tax_most_)[66,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[66,"Genus"], "__")[[1]][2]
    784. tax_table(ps.ng.tax_most_)[66,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[66,"Species"], "__")[[1]][2]
    785. tax_table(ps.ng.tax_most_)[67,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[67,"Domain"], "__")[[1]][2]
    786. tax_table(ps.ng.tax_most_)[67,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[67,"Phylum"], "__")[[1]][2]
    787. tax_table(ps.ng.tax_most_)[67,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[67,"Class"], "__")[[1]][2]
    788. tax_table(ps.ng.tax_most_)[67,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[67,"Order"], "__")[[1]][2]
    789. tax_table(ps.ng.tax_most_)[67,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[67,"Family"], "__")[[1]][2]
    790. tax_table(ps.ng.tax_most_)[67,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[67,"Genus"], "__")[[1]][2]
    791. tax_table(ps.ng.tax_most_)[67,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[67,"Species"], "__")[[1]][2]
    792. tax_table(ps.ng.tax_most_)[68,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[68,"Domain"], "__")[[1]][2]
    793. tax_table(ps.ng.tax_most_)[68,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[68,"Phylum"], "__")[[1]][2]
    794. tax_table(ps.ng.tax_most_)[68,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[68,"Class"], "__")[[1]][2]
    795. tax_table(ps.ng.tax_most_)[68,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[68,"Order"], "__")[[1]][2]
    796. tax_table(ps.ng.tax_most_)[68,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[68,"Family"], "__")[[1]][2]
    797. tax_table(ps.ng.tax_most_)[68,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[68,"Genus"], "__")[[1]][2]
    798. tax_table(ps.ng.tax_most_)[68,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[68,"Species"], "__")[[1]][2]
    799. tax_table(ps.ng.tax_most_)[69,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[69,"Domain"], "__")[[1]][2]
    800. tax_table(ps.ng.tax_most_)[69,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[69,"Phylum"], "__")[[1]][2]
    801. tax_table(ps.ng.tax_most_)[69,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[69,"Class"], "__")[[1]][2]
    802. tax_table(ps.ng.tax_most_)[69,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[69,"Order"], "__")[[1]][2]
    803. tax_table(ps.ng.tax_most_)[69,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[69,"Family"], "__")[[1]][2]
    804. tax_table(ps.ng.tax_most_)[69,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[69,"Genus"], "__")[[1]][2]
    805. tax_table(ps.ng.tax_most_)[69,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[69,"Species"], "__")[[1]][2]
    806. tax_table(ps.ng.tax_most_)[70,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[70,"Domain"], "__")[[1]][2]
    807. tax_table(ps.ng.tax_most_)[70,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[70,"Phylum"], "__")[[1]][2]
    808. tax_table(ps.ng.tax_most_)[70,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[70,"Class"], "__")[[1]][2]
    809. tax_table(ps.ng.tax_most_)[70,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[70,"Order"], "__")[[1]][2]
    810. tax_table(ps.ng.tax_most_)[70,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[70,"Family"], "__")[[1]][2]
    811. tax_table(ps.ng.tax_most_)[70,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[70,"Genus"], "__")[[1]][2]
    812. tax_table(ps.ng.tax_most_)[70,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[70,"Species"], "__")[[1]][2]
    813. tax_table(ps.ng.tax_most_)[71,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[71,"Domain"], "__")[[1]][2]
    814. tax_table(ps.ng.tax_most_)[71,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[71,"Phylum"], "__")[[1]][2]
    815. tax_table(ps.ng.tax_most_)[71,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[71,"Class"], "__")[[1]][2]
    816. tax_table(ps.ng.tax_most_)[71,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[71,"Order"], "__")[[1]][2]
    817. tax_table(ps.ng.tax_most_)[71,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[71,"Family"], "__")[[1]][2]
    818. tax_table(ps.ng.tax_most_)[71,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[71,"Genus"], "__")[[1]][2]
    819. tax_table(ps.ng.tax_most_)[71,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[71,"Species"], "__")[[1]][2]
    820. tax_table(ps.ng.tax_most_)[72,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[72,"Domain"], "__")[[1]][2]
    821. tax_table(ps.ng.tax_most_)[72,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[72,"Phylum"], "__")[[1]][2]
    822. tax_table(ps.ng.tax_most_)[72,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[72,"Class"], "__")[[1]][2]
    823. tax_table(ps.ng.tax_most_)[72,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[72,"Order"], "__")[[1]][2]
    824. tax_table(ps.ng.tax_most_)[72,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[72,"Family"], "__")[[1]][2]
    825. tax_table(ps.ng.tax_most_)[72,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[72,"Genus"], "__")[[1]][2]
    826. tax_table(ps.ng.tax_most_)[72,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[72,"Species"], "__")[[1]][2]
    827. tax_table(ps.ng.tax_most_)[73,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[73,"Domain"], "__")[[1]][2]
    828. tax_table(ps.ng.tax_most_)[73,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[73,"Phylum"], "__")[[1]][2]
    829. tax_table(ps.ng.tax_most_)[73,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[73,"Class"], "__")[[1]][2]
    830. tax_table(ps.ng.tax_most_)[73,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[73,"Order"], "__")[[1]][2]
    831. tax_table(ps.ng.tax_most_)[73,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[73,"Family"], "__")[[1]][2]
    832. tax_table(ps.ng.tax_most_)[73,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[73,"Genus"], "__")[[1]][2]
    833. tax_table(ps.ng.tax_most_)[73,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[73,"Species"], "__")[[1]][2]
    834. tax_table(ps.ng.tax_most_)[74,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[74,"Domain"], "__")[[1]][2]
    835. tax_table(ps.ng.tax_most_)[74,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[74,"Phylum"], "__")[[1]][2]
    836. tax_table(ps.ng.tax_most_)[74,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[74,"Class"], "__")[[1]][2]
    837. tax_table(ps.ng.tax_most_)[74,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[74,"Order"], "__")[[1]][2]
    838. tax_table(ps.ng.tax_most_)[74,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[74,"Family"], "__")[[1]][2]
    839. tax_table(ps.ng.tax_most_)[74,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[74,"Genus"], "__")[[1]][2]
    840. tax_table(ps.ng.tax_most_)[74,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[74,"Species"], "__")[[1]][2]
    841. tax_table(ps.ng.tax_most_)[75,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[75,"Domain"], "__")[[1]][2]
    842. tax_table(ps.ng.tax_most_)[75,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[75,"Phylum"], "__")[[1]][2]
    843. tax_table(ps.ng.tax_most_)[75,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[75,"Class"], "__")[[1]][2]
    844. tax_table(ps.ng.tax_most_)[75,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[75,"Order"], "__")[[1]][2]
    845. tax_table(ps.ng.tax_most_)[75,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[75,"Family"], "__")[[1]][2]
    846. tax_table(ps.ng.tax_most_)[75,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[75,"Genus"], "__")[[1]][2]
    847. tax_table(ps.ng.tax_most_)[75,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[75,"Species"], "__")[[1]][2]
    848. tax_table(ps.ng.tax_most_)[76,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[76,"Domain"], "__")[[1]][2]
    849. tax_table(ps.ng.tax_most_)[76,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[76,"Phylum"], "__")[[1]][2]
    850. tax_table(ps.ng.tax_most_)[76,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[76,"Class"], "__")[[1]][2]
    851. tax_table(ps.ng.tax_most_)[76,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[76,"Order"], "__")[[1]][2]
    852. tax_table(ps.ng.tax_most_)[76,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[76,"Family"], "__")[[1]][2]
    853. tax_table(ps.ng.tax_most_)[76,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[76,"Genus"], "__")[[1]][2]
    854. tax_table(ps.ng.tax_most_)[76,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[76,"Species"], "__")[[1]][2]
    855. tax_table(ps.ng.tax_most_)[77,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[77,"Domain"], "__")[[1]][2]
    856. tax_table(ps.ng.tax_most_)[77,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[77,"Phylum"], "__")[[1]][2]
    857. tax_table(ps.ng.tax_most_)[77,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[77,"Class"], "__")[[1]][2]
    858. tax_table(ps.ng.tax_most_)[77,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[77,"Order"], "__")[[1]][2]
    859. tax_table(ps.ng.tax_most_)[77,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[77,"Family"], "__")[[1]][2]
    860. tax_table(ps.ng.tax_most_)[77,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[77,"Genus"], "__")[[1]][2]
    861. tax_table(ps.ng.tax_most_)[77,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[77,"Species"], "__")[[1]][2]
    862. tax_table(ps.ng.tax_most_)[78,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[78,"Domain"], "__")[[1]][2]
    863. tax_table(ps.ng.tax_most_)[78,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[78,"Phylum"], "__")[[1]][2]
    864. tax_table(ps.ng.tax_most_)[78,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[78,"Class"], "__")[[1]][2]
    865. tax_table(ps.ng.tax_most_)[78,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[78,"Order"], "__")[[1]][2]
    866. tax_table(ps.ng.tax_most_)[78,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[78,"Family"], "__")[[1]][2]
    867. tax_table(ps.ng.tax_most_)[78,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[78,"Genus"], "__")[[1]][2]
    868. tax_table(ps.ng.tax_most_)[78,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[78,"Species"], "__")[[1]][2]
    869. tax_table(ps.ng.tax_most_)[79,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[79,"Domain"], "__")[[1]][2]
    870. tax_table(ps.ng.tax_most_)[79,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[79,"Phylum"], "__")[[1]][2]
    871. tax_table(ps.ng.tax_most_)[79,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[79,"Class"], "__")[[1]][2]
    872. tax_table(ps.ng.tax_most_)[79,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[79,"Order"], "__")[[1]][2]
    873. tax_table(ps.ng.tax_most_)[79,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[79,"Family"], "__")[[1]][2]
    874. tax_table(ps.ng.tax_most_)[79,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[79,"Genus"], "__")[[1]][2]
    875. tax_table(ps.ng.tax_most_)[79,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[79,"Species"], "__")[[1]][2]
    876. tax_table(ps.ng.tax_most_)[80,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[80,"Domain"], "__")[[1]][2]
    877. tax_table(ps.ng.tax_most_)[80,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[80,"Phylum"], "__")[[1]][2]
    878. tax_table(ps.ng.tax_most_)[80,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[80,"Class"], "__")[[1]][2]
    879. tax_table(ps.ng.tax_most_)[80,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[80,"Order"], "__")[[1]][2]
    880. tax_table(ps.ng.tax_most_)[80,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[80,"Family"], "__")[[1]][2]
    881. tax_table(ps.ng.tax_most_)[80,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[80,"Genus"], "__")[[1]][2]
    882. tax_table(ps.ng.tax_most_)[80,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[80,"Species"], "__")[[1]][2]
    883. tax_table(ps.ng.tax_most_)[81,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[81,"Domain"], "__")[[1]][2]
    884. tax_table(ps.ng.tax_most_)[81,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[81,"Phylum"], "__")[[1]][2]
    885. tax_table(ps.ng.tax_most_)[81,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[81,"Class"], "__")[[1]][2]
    886. tax_table(ps.ng.tax_most_)[81,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[81,"Order"], "__")[[1]][2]
    887. tax_table(ps.ng.tax_most_)[81,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[81,"Family"], "__")[[1]][2]
    888. tax_table(ps.ng.tax_most_)[81,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[81,"Genus"], "__")[[1]][2]
    889. tax_table(ps.ng.tax_most_)[81,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[81,"Species"], "__")[[1]][2]
    890. tax_table(ps.ng.tax_most_)[82,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[82,"Domain"], "__")[[1]][2]
    891. tax_table(ps.ng.tax_most_)[82,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[82,"Phylum"], "__")[[1]][2]
    892. tax_table(ps.ng.tax_most_)[82,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[82,"Class"], "__")[[1]][2]
    893. tax_table(ps.ng.tax_most_)[82,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[82,"Order"], "__")[[1]][2]
    894. tax_table(ps.ng.tax_most_)[82,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[82,"Family"], "__")[[1]][2]
    895. tax_table(ps.ng.tax_most_)[82,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[82,"Genus"], "__")[[1]][2]
    896. tax_table(ps.ng.tax_most_)[82,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[82,"Species"], "__")[[1]][2]
    897. tax_table(ps.ng.tax_most_)[83,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[83,"Domain"], "__")[[1]][2]
    898. tax_table(ps.ng.tax_most_)[83,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[83,"Phylum"], "__")[[1]][2]
    899. tax_table(ps.ng.tax_most_)[83,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[83,"Class"], "__")[[1]][2]
    900. tax_table(ps.ng.tax_most_)[83,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[83,"Order"], "__")[[1]][2]
    901. tax_table(ps.ng.tax_most_)[83,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[83,"Family"], "__")[[1]][2]
    902. tax_table(ps.ng.tax_most_)[83,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[83,"Genus"], "__")[[1]][2]
    903. tax_table(ps.ng.tax_most_)[83,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[83,"Species"], "__")[[1]][2]
    904. tax_table(ps.ng.tax_most_)[84,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[84,"Domain"], "__")[[1]][2]
    905. tax_table(ps.ng.tax_most_)[84,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[84,"Phylum"], "__")[[1]][2]
    906. tax_table(ps.ng.tax_most_)[84,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[84,"Class"], "__")[[1]][2]
    907. tax_table(ps.ng.tax_most_)[84,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[84,"Order"], "__")[[1]][2]
    908. tax_table(ps.ng.tax_most_)[84,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[84,"Family"], "__")[[1]][2]
    909. tax_table(ps.ng.tax_most_)[84,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[84,"Genus"], "__")[[1]][2]
    910. tax_table(ps.ng.tax_most_)[84,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[84,"Species"], "__")[[1]][2]
    911. tax_table(ps.ng.tax_most_)[85,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[85,"Domain"], "__")[[1]][2]
    912. tax_table(ps.ng.tax_most_)[85,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[85,"Phylum"], "__")[[1]][2]
    913. tax_table(ps.ng.tax_most_)[85,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[85,"Class"], "__")[[1]][2]
    914. tax_table(ps.ng.tax_most_)[85,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[85,"Order"], "__")[[1]][2]
    915. tax_table(ps.ng.tax_most_)[85,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[85,"Family"], "__")[[1]][2]
    916. tax_table(ps.ng.tax_most_)[85,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[85,"Genus"], "__")[[1]][2]
    917. tax_table(ps.ng.tax_most_)[85,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[85,"Species"], "__")[[1]][2]
    918. tax_table(ps.ng.tax_most_)[86,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[86,"Domain"], "__")[[1]][2]
    919. tax_table(ps.ng.tax_most_)[86,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[86,"Phylum"], "__")[[1]][2]
    920. tax_table(ps.ng.tax_most_)[86,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[86,"Class"], "__")[[1]][2]
    921. tax_table(ps.ng.tax_most_)[86,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[86,"Order"], "__")[[1]][2]
    922. tax_table(ps.ng.tax_most_)[86,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[86,"Family"], "__")[[1]][2]
    923. tax_table(ps.ng.tax_most_)[86,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[86,"Genus"], "__")[[1]][2]
    924. tax_table(ps.ng.tax_most_)[86,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[86,"Species"], "__")[[1]][2]
    925. tax_table(ps.ng.tax_most_)[87,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[87,"Domain"], "__")[[1]][2]
    926. tax_table(ps.ng.tax_most_)[87,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[87,"Phylum"], "__")[[1]][2]
    927. tax_table(ps.ng.tax_most_)[87,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[87,"Class"], "__")[[1]][2]
    928. tax_table(ps.ng.tax_most_)[87,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[87,"Order"], "__")[[1]][2]
    929. tax_table(ps.ng.tax_most_)[87,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[87,"Family"], "__")[[1]][2]
    930. tax_table(ps.ng.tax_most_)[87,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[87,"Genus"], "__")[[1]][2]
    931. tax_table(ps.ng.tax_most_)[87,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[87,"Species"], "__")[[1]][2]
    932. tax_table(ps.ng.tax_most_)[88,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[88,"Domain"], "__")[[1]][2]
    933. tax_table(ps.ng.tax_most_)[88,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[88,"Phylum"], "__")[[1]][2]
    934. tax_table(ps.ng.tax_most_)[88,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[88,"Class"], "__")[[1]][2]
    935. tax_table(ps.ng.tax_most_)[88,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[88,"Order"], "__")[[1]][2]
    936. tax_table(ps.ng.tax_most_)[88,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[88,"Family"], "__")[[1]][2]
    937. tax_table(ps.ng.tax_most_)[88,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[88,"Genus"], "__")[[1]][2]
    938. tax_table(ps.ng.tax_most_)[88,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[88,"Species"], "__")[[1]][2]
    939. tax_table(ps.ng.tax_most_)[89,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[89,"Domain"], "__")[[1]][2]
    940. tax_table(ps.ng.tax_most_)[89,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[89,"Phylum"], "__")[[1]][2]
    941. tax_table(ps.ng.tax_most_)[89,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[89,"Class"], "__")[[1]][2]
    942. tax_table(ps.ng.tax_most_)[89,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[89,"Order"], "__")[[1]][2]
    943. tax_table(ps.ng.tax_most_)[89,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[89,"Family"], "__")[[1]][2]
    944. tax_table(ps.ng.tax_most_)[89,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[89,"Genus"], "__")[[1]][2]
    945. tax_table(ps.ng.tax_most_)[89,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[89,"Species"], "__")[[1]][2]
    946. tax_table(ps.ng.tax_most_)[90,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[90,"Domain"], "__")[[1]][2]
    947. tax_table(ps.ng.tax_most_)[90,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[90,"Phylum"], "__")[[1]][2]
    948. tax_table(ps.ng.tax_most_)[90,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[90,"Class"], "__")[[1]][2]
    949. tax_table(ps.ng.tax_most_)[90,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[90,"Order"], "__")[[1]][2]
    950. tax_table(ps.ng.tax_most_)[90,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[90,"Family"], "__")[[1]][2]
    951. tax_table(ps.ng.tax_most_)[90,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[90,"Genus"], "__")[[1]][2]
    952. tax_table(ps.ng.tax_most_)[90,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[90,"Species"], "__")[[1]][2]
    953. tax_table(ps.ng.tax_most_)[91,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[91,"Domain"], "__")[[1]][2]
    954. tax_table(ps.ng.tax_most_)[91,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[91,"Phylum"], "__")[[1]][2]
    955. tax_table(ps.ng.tax_most_)[91,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[91,"Class"], "__")[[1]][2]
    956. tax_table(ps.ng.tax_most_)[91,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[91,"Order"], "__")[[1]][2]
    957. tax_table(ps.ng.tax_most_)[91,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[91,"Family"], "__")[[1]][2]
    958. tax_table(ps.ng.tax_most_)[91,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[91,"Genus"], "__")[[1]][2]
    959. tax_table(ps.ng.tax_most_)[91,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[91,"Species"], "__")[[1]][2]
    960. tax_table(ps.ng.tax_most_)[92,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[92,"Domain"], "__")[[1]][2]
    961. tax_table(ps.ng.tax_most_)[92,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[92,"Phylum"], "__")[[1]][2]
    962. tax_table(ps.ng.tax_most_)[92,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[92,"Class"], "__")[[1]][2]
    963. tax_table(ps.ng.tax_most_)[92,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[92,"Order"], "__")[[1]][2]
    964. tax_table(ps.ng.tax_most_)[92,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[92,"Family"], "__")[[1]][2]
    965. tax_table(ps.ng.tax_most_)[92,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[92,"Genus"], "__")[[1]][2]
    966. tax_table(ps.ng.tax_most_)[92,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[92,"Species"], "__")[[1]][2]
    967. tax_table(ps.ng.tax_most_)[93,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[93,"Domain"], "__")[[1]][2]
    968. tax_table(ps.ng.tax_most_)[93,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[93,"Phylum"], "__")[[1]][2]
    969. tax_table(ps.ng.tax_most_)[93,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[93,"Class"], "__")[[1]][2]
    970. tax_table(ps.ng.tax_most_)[93,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[93,"Order"], "__")[[1]][2]
    971. tax_table(ps.ng.tax_most_)[93,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[93,"Family"], "__")[[1]][2]
    972. tax_table(ps.ng.tax_most_)[93,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[93,"Genus"], "__")[[1]][2]
    973. tax_table(ps.ng.tax_most_)[93,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[93,"Species"], "__")[[1]][2]
    974. tax_table(ps.ng.tax_most_)[94,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[94,"Domain"], "__")[[1]][2]
    975. tax_table(ps.ng.tax_most_)[94,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[94,"Phylum"], "__")[[1]][2]
    976. tax_table(ps.ng.tax_most_)[94,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[94,"Class"], "__")[[1]][2]
    977. tax_table(ps.ng.tax_most_)[94,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[94,"Order"], "__")[[1]][2]
    978. tax_table(ps.ng.tax_most_)[94,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[94,"Family"], "__")[[1]][2]
    979. tax_table(ps.ng.tax_most_)[94,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[94,"Genus"], "__")[[1]][2]
    980. tax_table(ps.ng.tax_most_)[94,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[94,"Species"], "__")[[1]][2]
    981. tax_table(ps.ng.tax_most_)[95,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[95,"Domain"], "__")[[1]][2]
    982. tax_table(ps.ng.tax_most_)[95,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[95,"Phylum"], "__")[[1]][2]
    983. tax_table(ps.ng.tax_most_)[95,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[95,"Class"], "__")[[1]][2]
    984. tax_table(ps.ng.tax_most_)[95,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[95,"Order"], "__")[[1]][2]
    985. tax_table(ps.ng.tax_most_)[95,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[95,"Family"], "__")[[1]][2]
    986. tax_table(ps.ng.tax_most_)[95,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[95,"Genus"], "__")[[1]][2]
    987. tax_table(ps.ng.tax_most_)[95,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[95,"Species"], "__")[[1]][2]
    988. tax_table(ps.ng.tax_most_)[96,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[96,"Domain"], "__")[[1]][2]
    989. tax_table(ps.ng.tax_most_)[96,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[96,"Phylum"], "__")[[1]][2]
    990. tax_table(ps.ng.tax_most_)[96,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[96,"Class"], "__")[[1]][2]
    991. tax_table(ps.ng.tax_most_)[96,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[96,"Order"], "__")[[1]][2]
    992. tax_table(ps.ng.tax_most_)[96,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[96,"Family"], "__")[[1]][2]
    993. tax_table(ps.ng.tax_most_)[96,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[96,"Genus"], "__")[[1]][2]
    994. tax_table(ps.ng.tax_most_)[96,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[96,"Species"], "__")[[1]][2]
    995. tax_table(ps.ng.tax_most_)[97,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[97,"Domain"], "__")[[1]][2]
    996. tax_table(ps.ng.tax_most_)[97,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[97,"Phylum"], "__")[[1]][2]
    997. tax_table(ps.ng.tax_most_)[97,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[97,"Class"], "__")[[1]][2]
    998. tax_table(ps.ng.tax_most_)[97,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[97,"Order"], "__")[[1]][2]
    999. tax_table(ps.ng.tax_most_)[97,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[97,"Family"], "__")[[1]][2]
    1000. tax_table(ps.ng.tax_most_)[97,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[97,"Genus"], "__")[[1]][2]
    1001. tax_table(ps.ng.tax_most_)[97,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[97,"Species"], "__")[[1]][2]
    1002. tax_table(ps.ng.tax_most_)[98,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[98,"Domain"], "__")[[1]][2]
    1003. tax_table(ps.ng.tax_most_)[98,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[98,"Phylum"], "__")[[1]][2]
    1004. tax_table(ps.ng.tax_most_)[98,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[98,"Class"], "__")[[1]][2]
    1005. tax_table(ps.ng.tax_most_)[98,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[98,"Order"], "__")[[1]][2]
    1006. tax_table(ps.ng.tax_most_)[98,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[98,"Family"], "__")[[1]][2]
    1007. tax_table(ps.ng.tax_most_)[98,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[98,"Genus"], "__")[[1]][2]
    1008. tax_table(ps.ng.tax_most_)[98,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[98,"Species"], "__")[[1]][2]
    1009. tax_table(ps.ng.tax_most_)[99,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[99,"Domain"], "__")[[1]][2]
    1010. tax_table(ps.ng.tax_most_)[99,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[99,"Phylum"], "__")[[1]][2]
    1011. tax_table(ps.ng.tax_most_)[99,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[99,"Class"], "__")[[1]][2]
    1012. tax_table(ps.ng.tax_most_)[99,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[99,"Order"], "__")[[1]][2]
    1013. tax_table(ps.ng.tax_most_)[99,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[99,"Family"], "__")[[1]][2]
    1014. tax_table(ps.ng.tax_most_)[99,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[99,"Genus"], "__")[[1]][2]
    1015. tax_table(ps.ng.tax_most_)[99,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[99,"Species"], "__")[[1]][2]
    1016. tax_table(ps.ng.tax_most_)[100,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[100,"Domain"], "__")[[1]][2]
    1017. tax_table(ps.ng.tax_most_)[100,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[100,"Phylum"], "__")[[1]][2]
    1018. tax_table(ps.ng.tax_most_)[100,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[100,"Class"], "__")[[1]][2]
    1019. tax_table(ps.ng.tax_most_)[100,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[100,"Order"], "__")[[1]][2]
    1020. tax_table(ps.ng.tax_most_)[100,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[100,"Family"], "__")[[1]][2]
    1021. tax_table(ps.ng.tax_most_)[100,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[100,"Genus"], "__")[[1]][2]
    1022. tax_table(ps.ng.tax_most_)[100,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[100,"Species"], "__")[[1]][2]
    1023. tax_table(ps.ng.tax_most_)[101,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[101,"Domain"], "__")[[1]][2]
    1024. tax_table(ps.ng.tax_most_)[101,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[101,"Phylum"], "__")[[1]][2]
    1025. tax_table(ps.ng.tax_most_)[101,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[101,"Class"], "__")[[1]][2]
    1026. tax_table(ps.ng.tax_most_)[101,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[101,"Order"], "__")[[1]][2]
    1027. tax_table(ps.ng.tax_most_)[101,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[101,"Family"], "__")[[1]][2]
    1028. tax_table(ps.ng.tax_most_)[101,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[101,"Genus"], "__")[[1]][2]
    1029. tax_table(ps.ng.tax_most_)[101,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[101,"Species"], "__")[[1]][2]
    1030. tax_table(ps.ng.tax_most_)[102,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[102,"Domain"], "__")[[1]][2]
    1031. tax_table(ps.ng.tax_most_)[102,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[102,"Phylum"], "__")[[1]][2]
    1032. tax_table(ps.ng.tax_most_)[102,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[102,"Class"], "__")[[1]][2]
    1033. tax_table(ps.ng.tax_most_)[102,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[102,"Order"], "__")[[1]][2]
    1034. tax_table(ps.ng.tax_most_)[102,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[102,"Family"], "__")[[1]][2]
    1035. tax_table(ps.ng.tax_most_)[102,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[102,"Genus"], "__")[[1]][2]
    1036. tax_table(ps.ng.tax_most_)[102,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[102,"Species"], "__")[[1]][2]
    1037. tax_table(ps.ng.tax_most_)[103,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[103,"Domain"], "__")[[1]][2]
    1038. tax_table(ps.ng.tax_most_)[103,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[103,"Phylum"], "__")[[1]][2]
    1039. tax_table(ps.ng.tax_most_)[103,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[103,"Class"], "__")[[1]][2]
    1040. tax_table(ps.ng.tax_most_)[103,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[103,"Order"], "__")[[1]][2]
    1041. tax_table(ps.ng.tax_most_)[103,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[103,"Family"], "__")[[1]][2]
    1042. tax_table(ps.ng.tax_most_)[103,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[103,"Genus"], "__")[[1]][2]
    1043. tax_table(ps.ng.tax_most_)[103,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[103,"Species"], "__")[[1]][2]
    1044. tax_table(ps.ng.tax_most_)[104,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[104,"Domain"], "__")[[1]][2]
    1045. tax_table(ps.ng.tax_most_)[104,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[104,"Phylum"], "__")[[1]][2]
    1046. tax_table(ps.ng.tax_most_)[104,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[104,"Class"], "__")[[1]][2]
    1047. tax_table(ps.ng.tax_most_)[104,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[104,"Order"], "__")[[1]][2]
    1048. tax_table(ps.ng.tax_most_)[104,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[104,"Family"], "__")[[1]][2]
    1049. tax_table(ps.ng.tax_most_)[104,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[104,"Genus"], "__")[[1]][2]
    1050. tax_table(ps.ng.tax_most_)[104,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[104,"Species"], "__")[[1]][2]
    1051. tax_table(ps.ng.tax_most_)[105,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[105,"Domain"], "__")[[1]][2]
    1052. tax_table(ps.ng.tax_most_)[105,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[105,"Phylum"], "__")[[1]][2]
    1053. tax_table(ps.ng.tax_most_)[105,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[105,"Class"], "__")[[1]][2]
    1054. tax_table(ps.ng.tax_most_)[105,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[105,"Order"], "__")[[1]][2]
    1055. tax_table(ps.ng.tax_most_)[105,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[105,"Family"], "__")[[1]][2]
    1056. tax_table(ps.ng.tax_most_)[105,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[105,"Genus"], "__")[[1]][2]
    1057. tax_table(ps.ng.tax_most_)[105,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[105,"Species"], "__")[[1]][2]
    1058. tax_table(ps.ng.tax_most_)[106,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[106,"Domain"], "__")[[1]][2]
    1059. tax_table(ps.ng.tax_most_)[106,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[106,"Phylum"], "__")[[1]][2]
    1060. tax_table(ps.ng.tax_most_)[106,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[106,"Class"], "__")[[1]][2]
    1061. tax_table(ps.ng.tax_most_)[106,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[106,"Order"], "__")[[1]][2]
    1062. tax_table(ps.ng.tax_most_)[106,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[106,"Family"], "__")[[1]][2]
    1063. tax_table(ps.ng.tax_most_)[106,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[106,"Genus"], "__")[[1]][2]
    1064. tax_table(ps.ng.tax_most_)[106,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[106,"Species"], "__")[[1]][2]
    1065. tax_table(ps.ng.tax_most_)[107,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[107,"Domain"], "__")[[1]][2]
    1066. tax_table(ps.ng.tax_most_)[107,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[107,"Phylum"], "__")[[1]][2]
    1067. tax_table(ps.ng.tax_most_)[107,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[107,"Class"], "__")[[1]][2]
    1068. tax_table(ps.ng.tax_most_)[107,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[107,"Order"], "__")[[1]][2]
    1069. tax_table(ps.ng.tax_most_)[107,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[107,"Family"], "__")[[1]][2]
    1070. tax_table(ps.ng.tax_most_)[107,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[107,"Genus"], "__")[[1]][2]
    1071. tax_table(ps.ng.tax_most_)[107,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[107,"Species"], "__")[[1]][2]
    1072. tax_table(ps.ng.tax_most_)[108,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[108,"Domain"], "__")[[1]][2]
    1073. tax_table(ps.ng.tax_most_)[108,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[108,"Phylum"], "__")[[1]][2]
    1074. tax_table(ps.ng.tax_most_)[108,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[108,"Class"], "__")[[1]][2]
    1075. tax_table(ps.ng.tax_most_)[108,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[108,"Order"], "__")[[1]][2]
    1076. tax_table(ps.ng.tax_most_)[108,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[108,"Family"], "__")[[1]][2]
    1077. tax_table(ps.ng.tax_most_)[108,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[108,"Genus"], "__")[[1]][2]
    1078. tax_table(ps.ng.tax_most_)[108,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[108,"Species"], "__")[[1]][2]
    1079. tax_table(ps.ng.tax_most_)[109,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[109,"Domain"], "__")[[1]][2]
    1080. tax_table(ps.ng.tax_most_)[109,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[109,"Phylum"], "__")[[1]][2]
    1081. tax_table(ps.ng.tax_most_)[109,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[109,"Class"], "__")[[1]][2]
    1082. tax_table(ps.ng.tax_most_)[109,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[109,"Order"], "__")[[1]][2]
    1083. tax_table(ps.ng.tax_most_)[109,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[109,"Family"], "__")[[1]][2]
    1084. tax_table(ps.ng.tax_most_)[109,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[109,"Genus"], "__")[[1]][2]
    1085. tax_table(ps.ng.tax_most_)[109,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[109,"Species"], "__")[[1]][2]
    1086. tax_table(ps.ng.tax_most_)[110,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[110,"Domain"], "__")[[1]][2]
    1087. tax_table(ps.ng.tax_most_)[110,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[110,"Phylum"], "__")[[1]][2]
    1088. tax_table(ps.ng.tax_most_)[110,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[110,"Class"], "__")[[1]][2]
    1089. tax_table(ps.ng.tax_most_)[110,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[110,"Order"], "__")[[1]][2]
    1090. tax_table(ps.ng.tax_most_)[110,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[110,"Family"], "__")[[1]][2]
    1091. tax_table(ps.ng.tax_most_)[110,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[110,"Genus"], "__")[[1]][2]
    1092. tax_table(ps.ng.tax_most_)[110,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[110,"Species"], "__")[[1]][2]
    1093. tax_table(ps.ng.tax_most_)[111,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[111,"Domain"], "__")[[1]][2]
    1094. tax_table(ps.ng.tax_most_)[111,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[111,"Phylum"], "__")[[1]][2]
    1095. tax_table(ps.ng.tax_most_)[111,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[111,"Class"], "__")[[1]][2]
    1096. tax_table(ps.ng.tax_most_)[111,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[111,"Order"], "__")[[1]][2]
    1097. tax_table(ps.ng.tax_most_)[111,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[111,"Family"], "__")[[1]][2]
    1098. tax_table(ps.ng.tax_most_)[111,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[111,"Genus"], "__")[[1]][2]
    1099. tax_table(ps.ng.tax_most_)[111,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[111,"Species"], "__")[[1]][2]
    1100. tax_table(ps.ng.tax_most_)[112,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[112,"Domain"], "__")[[1]][2]
    1101. tax_table(ps.ng.tax_most_)[112,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[112,"Phylum"], "__")[[1]][2]
    1102. tax_table(ps.ng.tax_most_)[112,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[112,"Class"], "__")[[1]][2]
    1103. tax_table(ps.ng.tax_most_)[112,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[112,"Order"], "__")[[1]][2]
    1104. tax_table(ps.ng.tax_most_)[112,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[112,"Family"], "__")[[1]][2]
    1105. tax_table(ps.ng.tax_most_)[112,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[112,"Genus"], "__")[[1]][2]
    1106. tax_table(ps.ng.tax_most_)[112,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[112,"Species"], "__")[[1]][2]
    1107. tax_table(ps.ng.tax_most_)[113,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[113,"Domain"], "__")[[1]][2]
    1108. tax_table(ps.ng.tax_most_)[113,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[113,"Phylum"], "__")[[1]][2]
    1109. tax_table(ps.ng.tax_most_)[113,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[113,"Class"], "__")[[1]][2]
    1110. tax_table(ps.ng.tax_most_)[113,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[113,"Order"], "__")[[1]][2]
    1111. tax_table(ps.ng.tax_most_)[113,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[113,"Family"], "__")[[1]][2]
    1112. tax_table(ps.ng.tax_most_)[113,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[113,"Genus"], "__")[[1]][2]
    1113. tax_table(ps.ng.tax_most_)[113,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[113,"Species"], "__")[[1]][2]
    1114. tax_table(ps.ng.tax_most_)[114,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[114,"Domain"], "__")[[1]][2]
    1115. tax_table(ps.ng.tax_most_)[114,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[114,"Phylum"], "__")[[1]][2]
    1116. tax_table(ps.ng.tax_most_)[114,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[114,"Class"], "__")[[1]][2]
    1117. tax_table(ps.ng.tax_most_)[114,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[114,"Order"], "__")[[1]][2]
    1118. tax_table(ps.ng.tax_most_)[114,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[114,"Family"], "__")[[1]][2]
    1119. tax_table(ps.ng.tax_most_)[114,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[114,"Genus"], "__")[[1]][2]
    1120. tax_table(ps.ng.tax_most_)[114,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[114,"Species"], "__")[[1]][2]
    1121. tax_table(ps.ng.tax_most_)[115,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[115,"Domain"], "__")[[1]][2]
    1122. tax_table(ps.ng.tax_most_)[115,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[115,"Phylum"], "__")[[1]][2]
    1123. tax_table(ps.ng.tax_most_)[115,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[115,"Class"], "__")[[1]][2]
    1124. tax_table(ps.ng.tax_most_)[115,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[115,"Order"], "__")[[1]][2]
    1125. tax_table(ps.ng.tax_most_)[115,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[115,"Family"], "__")[[1]][2]
    1126. tax_table(ps.ng.tax_most_)[115,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[115,"Genus"], "__")[[1]][2]
    1127. tax_table(ps.ng.tax_most_)[115,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[115,"Species"], "__")[[1]][2]
    1128. tax_table(ps.ng.tax_most_)[116,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[116,"Domain"], "__")[[1]][2]
    1129. tax_table(ps.ng.tax_most_)[116,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[116,"Phylum"], "__")[[1]][2]
    1130. tax_table(ps.ng.tax_most_)[116,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[116,"Class"], "__")[[1]][2]
    1131. tax_table(ps.ng.tax_most_)[116,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[116,"Order"], "__")[[1]][2]
    1132. tax_table(ps.ng.tax_most_)[116,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[116,"Family"], "__")[[1]][2]
    1133. tax_table(ps.ng.tax_most_)[116,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[116,"Genus"], "__")[[1]][2]
    1134. tax_table(ps.ng.tax_most_)[116,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[116,"Species"], "__")[[1]][2]
    1135. tax_table(ps.ng.tax_most_)[117,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[117,"Domain"], "__")[[1]][2]
    1136. tax_table(ps.ng.tax_most_)[117,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[117,"Phylum"], "__")[[1]][2]
    1137. tax_table(ps.ng.tax_most_)[117,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[117,"Class"], "__")[[1]][2]
    1138. tax_table(ps.ng.tax_most_)[117,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[117,"Order"], "__")[[1]][2]
    1139. tax_table(ps.ng.tax_most_)[117,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[117,"Family"], "__")[[1]][2]
    1140. tax_table(ps.ng.tax_most_)[117,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[117,"Genus"], "__")[[1]][2]
    1141. tax_table(ps.ng.tax_most_)[117,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[117,"Species"], "__")[[1]][2]
    1142. tax_table(ps.ng.tax_most_)[118,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[118,"Domain"], "__")[[1]][2]
    1143. tax_table(ps.ng.tax_most_)[118,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[118,"Phylum"], "__")[[1]][2]
    1144. tax_table(ps.ng.tax_most_)[118,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[118,"Class"], "__")[[1]][2]
    1145. tax_table(ps.ng.tax_most_)[118,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[118,"Order"], "__")[[1]][2]
    1146. tax_table(ps.ng.tax_most_)[118,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[118,"Family"], "__")[[1]][2]
    1147. tax_table(ps.ng.tax_most_)[118,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[118,"Genus"], "__")[[1]][2]
    1148. tax_table(ps.ng.tax_most_)[118,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[118,"Species"], "__")[[1]][2]
    1149. tax_table(ps.ng.tax_most_)[119,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[119,"Domain"], "__")[[1]][2]
    1150. tax_table(ps.ng.tax_most_)[119,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[119,"Phylum"], "__")[[1]][2]
    1151. tax_table(ps.ng.tax_most_)[119,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[119,"Class"], "__")[[1]][2]
    1152. tax_table(ps.ng.tax_most_)[119,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[119,"Order"], "__")[[1]][2]
    1153. tax_table(ps.ng.tax_most_)[119,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[119,"Family"], "__")[[1]][2]
    1154. tax_table(ps.ng.tax_most_)[119,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[119,"Genus"], "__")[[1]][2]
    1155. tax_table(ps.ng.tax_most_)[119,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[119,"Species"], "__")[[1]][2]
    1156. tax_table(ps.ng.tax_most_)[120,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[120,"Domain"], "__")[[1]][2]
    1157. tax_table(ps.ng.tax_most_)[120,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[120,"Phylum"], "__")[[1]][2]
    1158. tax_table(ps.ng.tax_most_)[120,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[120,"Class"], "__")[[1]][2]
    1159. tax_table(ps.ng.tax_most_)[120,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[120,"Order"], "__")[[1]][2]
    1160. tax_table(ps.ng.tax_most_)[120,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[120,"Family"], "__")[[1]][2]
    1161. tax_table(ps.ng.tax_most_)[120,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[120,"Genus"], "__")[[1]][2]
    1162. tax_table(ps.ng.tax_most_)[120,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[120,"Species"], "__")[[1]][2]
    1163. tax_table(ps.ng.tax_most_)[121,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[121,"Domain"], "__")[[1]][2]
    1164. tax_table(ps.ng.tax_most_)[121,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[121,"Phylum"], "__")[[1]][2]
    1165. tax_table(ps.ng.tax_most_)[121,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[121,"Class"], "__")[[1]][2]
    1166. tax_table(ps.ng.tax_most_)[121,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[121,"Order"], "__")[[1]][2]
    1167. tax_table(ps.ng.tax_most_)[121,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[121,"Family"], "__")[[1]][2]
    1168. tax_table(ps.ng.tax_most_)[121,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[121,"Genus"], "__")[[1]][2]
    1169. tax_table(ps.ng.tax_most_)[121,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[121,"Species"], "__")[[1]][2]
    1170. tax_table(ps.ng.tax_most_)[122,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[122,"Domain"], "__")[[1]][2]
    1171. tax_table(ps.ng.tax_most_)[122,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[122,"Phylum"], "__")[[1]][2]
    1172. tax_table(ps.ng.tax_most_)[122,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[122,"Class"], "__")[[1]][2]
    1173. tax_table(ps.ng.tax_most_)[122,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[122,"Order"], "__")[[1]][2]
    1174. tax_table(ps.ng.tax_most_)[122,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[122,"Family"], "__")[[1]][2]
    1175. tax_table(ps.ng.tax_most_)[122,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[122,"Genus"], "__")[[1]][2]
    1176. tax_table(ps.ng.tax_most_)[122,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[122,"Species"], "__")[[1]][2]
    1177. tax_table(ps.ng.tax_most_)[123,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[123,"Domain"], "__")[[1]][2]
    1178. tax_table(ps.ng.tax_most_)[123,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[123,"Phylum"], "__")[[1]][2]
    1179. tax_table(ps.ng.tax_most_)[123,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[123,"Class"], "__")[[1]][2]
    1180. tax_table(ps.ng.tax_most_)[123,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[123,"Order"], "__")[[1]][2]
    1181. tax_table(ps.ng.tax_most_)[123,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[123,"Family"], "__")[[1]][2]
    1182. tax_table(ps.ng.tax_most_)[123,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[123,"Genus"], "__")[[1]][2]
    1183. tax_table(ps.ng.tax_most_)[123,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[123,"Species"], "__")[[1]][2]
    1184. tax_table(ps.ng.tax_most_)[124,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[124,"Domain"], "__")[[1]][2]
    1185. tax_table(ps.ng.tax_most_)[124,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[124,"Phylum"], "__")[[1]][2]
    1186. tax_table(ps.ng.tax_most_)[124,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[124,"Class"], "__")[[1]][2]
    1187. tax_table(ps.ng.tax_most_)[124,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[124,"Order"], "__")[[1]][2]
    1188. tax_table(ps.ng.tax_most_)[124,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[124,"Family"], "__")[[1]][2]
    1189. tax_table(ps.ng.tax_most_)[124,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[124,"Genus"], "__")[[1]][2]
    1190. tax_table(ps.ng.tax_most_)[124,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[124,"Species"], "__")[[1]][2]
    1191. tax_table(ps.ng.tax_most_)[125,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[125,"Domain"], "__")[[1]][2]
    1192. tax_table(ps.ng.tax_most_)[125,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[125,"Phylum"], "__")[[1]][2]
    1193. tax_table(ps.ng.tax_most_)[125,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[125,"Class"], "__")[[1]][2]
    1194. tax_table(ps.ng.tax_most_)[125,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[125,"Order"], "__")[[1]][2]
    1195. tax_table(ps.ng.tax_most_)[125,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[125,"Family"], "__")[[1]][2]
    1196. tax_table(ps.ng.tax_most_)[125,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[125,"Genus"], "__")[[1]][2]
    1197. tax_table(ps.ng.tax_most_)[125,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[125,"Species"], "__")[[1]][2]
    1198. tax_table(ps.ng.tax_most_)[126,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[126,"Domain"], "__")[[1]][2]
    1199. tax_table(ps.ng.tax_most_)[126,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[126,"Phylum"], "__")[[1]][2]
    1200. tax_table(ps.ng.tax_most_)[126,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[126,"Class"], "__")[[1]][2]
    1201. tax_table(ps.ng.tax_most_)[126,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[126,"Order"], "__")[[1]][2]
    1202. tax_table(ps.ng.tax_most_)[126,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[126,"Family"], "__")[[1]][2]
    1203. tax_table(ps.ng.tax_most_)[126,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[126,"Genus"], "__")[[1]][2]
    1204. tax_table(ps.ng.tax_most_)[126,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[126,"Species"], "__")[[1]][2]
    1205. tax_table(ps.ng.tax_most_)[127,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[127,"Domain"], "__")[[1]][2]
    1206. tax_table(ps.ng.tax_most_)[127,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[127,"Phylum"], "__")[[1]][2]
    1207. tax_table(ps.ng.tax_most_)[127,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[127,"Class"], "__")[[1]][2]
    1208. tax_table(ps.ng.tax_most_)[127,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[127,"Order"], "__")[[1]][2]
    1209. tax_table(ps.ng.tax_most_)[127,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[127,"Family"], "__")[[1]][2]
    1210. tax_table(ps.ng.tax_most_)[127,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[127,"Genus"], "__")[[1]][2]
    1211. tax_table(ps.ng.tax_most_)[127,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[127,"Species"], "__")[[1]][2]
    1212. tax_table(ps.ng.tax_most_)[128,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[128,"Domain"], "__")[[1]][2]
    1213. tax_table(ps.ng.tax_most_)[128,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[128,"Phylum"], "__")[[1]][2]
    1214. tax_table(ps.ng.tax_most_)[128,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[128,"Class"], "__")[[1]][2]
    1215. tax_table(ps.ng.tax_most_)[128,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[128,"Order"], "__")[[1]][2]
    1216. tax_table(ps.ng.tax_most_)[128,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[128,"Family"], "__")[[1]][2]
    1217. tax_table(ps.ng.tax_most_)[128,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[128,"Genus"], "__")[[1]][2]
    1218. tax_table(ps.ng.tax_most_)[128,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[128,"Species"], "__")[[1]][2]
    1219. tax_table(ps.ng.tax_most_)[129,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[129,"Domain"], "__")[[1]][2]
    1220. tax_table(ps.ng.tax_most_)[129,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[129,"Phylum"], "__")[[1]][2]
    1221. tax_table(ps.ng.tax_most_)[129,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[129,"Class"], "__")[[1]][2]
    1222. tax_table(ps.ng.tax_most_)[129,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[129,"Order"], "__")[[1]][2]
    1223. tax_table(ps.ng.tax_most_)[129,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[129,"Family"], "__")[[1]][2]
    1224. tax_table(ps.ng.tax_most_)[129,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[129,"Genus"], "__")[[1]][2]
    1225. tax_table(ps.ng.tax_most_)[129,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[129,"Species"], "__")[[1]][2]
    1226. tax_table(ps.ng.tax_most_)[130,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[130,"Domain"], "__")[[1]][2]
    1227. tax_table(ps.ng.tax_most_)[130,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[130,"Phylum"], "__")[[1]][2]
    1228. tax_table(ps.ng.tax_most_)[130,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[130,"Class"], "__")[[1]][2]
    1229. tax_table(ps.ng.tax_most_)[130,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[130,"Order"], "__")[[1]][2]
    1230. tax_table(ps.ng.tax_most_)[130,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[130,"Family"], "__")[[1]][2]
    1231. tax_table(ps.ng.tax_most_)[130,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[130,"Genus"], "__")[[1]][2]
    1232. tax_table(ps.ng.tax_most_)[130,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[130,"Species"], "__")[[1]][2]
    1233. tax_table(ps.ng.tax_most_)[131,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[131,"Domain"], "__")[[1]][2]
    1234. tax_table(ps.ng.tax_most_)[131,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[131,"Phylum"], "__")[[1]][2]
    1235. tax_table(ps.ng.tax_most_)[131,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[131,"Class"], "__")[[1]][2]
    1236. tax_table(ps.ng.tax_most_)[131,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[131,"Order"], "__")[[1]][2]
    1237. tax_table(ps.ng.tax_most_)[131,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[131,"Family"], "__")[[1]][2]
    1238. tax_table(ps.ng.tax_most_)[131,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[131,"Genus"], "__")[[1]][2]
    1239. tax_table(ps.ng.tax_most_)[131,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[131,"Species"], "__")[[1]][2]
    1240. tax_table(ps.ng.tax_most_)[132,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[132,"Domain"], "__")[[1]][2]
    1241. tax_table(ps.ng.tax_most_)[132,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[132,"Phylum"], "__")[[1]][2]
    1242. tax_table(ps.ng.tax_most_)[132,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[132,"Class"], "__")[[1]][2]
    1243. tax_table(ps.ng.tax_most_)[132,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[132,"Order"], "__")[[1]][2]
    1244. tax_table(ps.ng.tax_most_)[132,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[132,"Family"], "__")[[1]][2]
    1245. tax_table(ps.ng.tax_most_)[132,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[132,"Genus"], "__")[[1]][2]
    1246. tax_table(ps.ng.tax_most_)[132,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[132,"Species"], "__")[[1]][2]
    1247. tax_table(ps.ng.tax_most_)[133,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[133,"Domain"], "__")[[1]][2]
    1248. tax_table(ps.ng.tax_most_)[133,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[133,"Phylum"], "__")[[1]][2]
    1249. tax_table(ps.ng.tax_most_)[133,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[133,"Class"], "__")[[1]][2]
    1250. tax_table(ps.ng.tax_most_)[133,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[133,"Order"], "__")[[1]][2]
    1251. tax_table(ps.ng.tax_most_)[133,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[133,"Family"], "__")[[1]][2]
    1252. tax_table(ps.ng.tax_most_)[133,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[133,"Genus"], "__")[[1]][2]
    1253. tax_table(ps.ng.tax_most_)[133,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[133,"Species"], "__")[[1]][2]
    1254. tax_table(ps.ng.tax_most_)[134,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[134,"Domain"], "__")[[1]][2]
    1255. tax_table(ps.ng.tax_most_)[134,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[134,"Phylum"], "__")[[1]][2]
    1256. tax_table(ps.ng.tax_most_)[134,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[134,"Class"], "__")[[1]][2]
    1257. tax_table(ps.ng.tax_most_)[134,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[134,"Order"], "__")[[1]][2]
    1258. tax_table(ps.ng.tax_most_)[134,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[134,"Family"], "__")[[1]][2]
    1259. tax_table(ps.ng.tax_most_)[134,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[134,"Genus"], "__")[[1]][2]
    1260. tax_table(ps.ng.tax_most_)[134,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[134,"Species"], "__")[[1]][2]
    1261. tax_table(ps.ng.tax_most_)[135,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[135,"Domain"], "__")[[1]][2]
    1262. tax_table(ps.ng.tax_most_)[135,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[135,"Phylum"], "__")[[1]][2]
    1263. tax_table(ps.ng.tax_most_)[135,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[135,"Class"], "__")[[1]][2]
    1264. tax_table(ps.ng.tax_most_)[135,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[135,"Order"], "__")[[1]][2]
    1265. tax_table(ps.ng.tax_most_)[135,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[135,"Family"], "__")[[1]][2]
    1266. tax_table(ps.ng.tax_most_)[135,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[135,"Genus"], "__")[[1]][2]
    1267. tax_table(ps.ng.tax_most_)[135,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[135,"Species"], "__")[[1]][2]
    1268. tax_table(ps.ng.tax_most_)[136,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[136,"Domain"], "__")[[1]][2]
    1269. tax_table(ps.ng.tax_most_)[136,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[136,"Phylum"], "__")[[1]][2]
    1270. tax_table(ps.ng.tax_most_)[136,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[136,"Class"], "__")[[1]][2]
    1271. tax_table(ps.ng.tax_most_)[136,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[136,"Order"], "__")[[1]][2]
    1272. tax_table(ps.ng.tax_most_)[136,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[136,"Family"], "__")[[1]][2]
    1273. tax_table(ps.ng.tax_most_)[136,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[136,"Genus"], "__")[[1]][2]
    1274. tax_table(ps.ng.tax_most_)[136,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[136,"Species"], "__")[[1]][2]
    1275. tax_table(ps.ng.tax_most_)[137,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[137,"Domain"], "__")[[1]][2]
    1276. tax_table(ps.ng.tax_most_)[137,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[137,"Phylum"], "__")[[1]][2]
    1277. tax_table(ps.ng.tax_most_)[137,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[137,"Class"], "__")[[1]][2]
    1278. tax_table(ps.ng.tax_most_)[137,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[137,"Order"], "__")[[1]][2]
    1279. tax_table(ps.ng.tax_most_)[137,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[137,"Family"], "__")[[1]][2]
    1280. tax_table(ps.ng.tax_most_)[137,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[137,"Genus"], "__")[[1]][2]
    1281. tax_table(ps.ng.tax_most_)[137,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[137,"Species"], "__")[[1]][2]
    1282. tax_table(ps.ng.tax_most_)[138,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[138,"Domain"], "__")[[1]][2]
    1283. tax_table(ps.ng.tax_most_)[138,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[138,"Phylum"], "__")[[1]][2]
    1284. tax_table(ps.ng.tax_most_)[138,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[138,"Class"], "__")[[1]][2]
    1285. tax_table(ps.ng.tax_most_)[138,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[138,"Order"], "__")[[1]][2]
    1286. tax_table(ps.ng.tax_most_)[138,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[138,"Family"], "__")[[1]][2]
    1287. tax_table(ps.ng.tax_most_)[138,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[138,"Genus"], "__")[[1]][2]
    1288. tax_table(ps.ng.tax_most_)[138,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[138,"Species"], "__")[[1]][2]
    1289. tax_table(ps.ng.tax_most_)[139,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[139,"Domain"], "__")[[1]][2]
    1290. tax_table(ps.ng.tax_most_)[139,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[139,"Phylum"], "__")[[1]][2]
    1291. tax_table(ps.ng.tax_most_)[139,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[139,"Class"], "__")[[1]][2]
    1292. tax_table(ps.ng.tax_most_)[139,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[139,"Order"], "__")[[1]][2]
    1293. tax_table(ps.ng.tax_most_)[139,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[139,"Family"], "__")[[1]][2]
    1294. tax_table(ps.ng.tax_most_)[139,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[139,"Genus"], "__")[[1]][2]
    1295. tax_table(ps.ng.tax_most_)[139,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[139,"Species"], "__")[[1]][2]
    1296. tax_table(ps.ng.tax_most_)[140,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[140,"Domain"], "__")[[1]][2]
    1297. tax_table(ps.ng.tax_most_)[140,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[140,"Phylum"], "__")[[1]][2]
    1298. tax_table(ps.ng.tax_most_)[140,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[140,"Class"], "__")[[1]][2]
    1299. tax_table(ps.ng.tax_most_)[140,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[140,"Order"], "__")[[1]][2]
    1300. tax_table(ps.ng.tax_most_)[140,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[140,"Family"], "__")[[1]][2]
    1301. tax_table(ps.ng.tax_most_)[140,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[140,"Genus"], "__")[[1]][2]
    1302. tax_table(ps.ng.tax_most_)[140,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[140,"Species"], "__")[[1]][2]
    1303. tax_table(ps.ng.tax_most_)[141,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[141,"Domain"], "__")[[1]][2]
    1304. tax_table(ps.ng.tax_most_)[141,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[141,"Phylum"], "__")[[1]][2]
    1305. tax_table(ps.ng.tax_most_)[141,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[141,"Class"], "__")[[1]][2]
    1306. tax_table(ps.ng.tax_most_)[141,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[141,"Order"], "__")[[1]][2]
    1307. tax_table(ps.ng.tax_most_)[141,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[141,"Family"], "__")[[1]][2]
    1308. tax_table(ps.ng.tax_most_)[141,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[141,"Genus"], "__")[[1]][2]
    1309. tax_table(ps.ng.tax_most_)[141,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[141,"Species"], "__")[[1]][2]
    1310. tax_table(ps.ng.tax_most_)[142,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[142,"Domain"], "__")[[1]][2]
    1311. tax_table(ps.ng.tax_most_)[142,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[142,"Phylum"], "__")[[1]][2]
    1312. tax_table(ps.ng.tax_most_)[142,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[142,"Class"], "__")[[1]][2]
    1313. tax_table(ps.ng.tax_most_)[142,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[142,"Order"], "__")[[1]][2]
    1314. tax_table(ps.ng.tax_most_)[142,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[142,"Family"], "__")[[1]][2]
    1315. tax_table(ps.ng.tax_most_)[142,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[142,"Genus"], "__")[[1]][2]
    1316. tax_table(ps.ng.tax_most_)[142,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[142,"Species"], "__")[[1]][2]
    1317. tax_table(ps.ng.tax_most_)[143,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[143,"Domain"], "__")[[1]][2]
    1318. tax_table(ps.ng.tax_most_)[143,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[143,"Phylum"], "__")[[1]][2]
    1319. tax_table(ps.ng.tax_most_)[143,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[143,"Class"], "__")[[1]][2]
    1320. tax_table(ps.ng.tax_most_)[143,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[143,"Order"], "__")[[1]][2]
    1321. tax_table(ps.ng.tax_most_)[143,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[143,"Family"], "__")[[1]][2]
    1322. tax_table(ps.ng.tax_most_)[143,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[143,"Genus"], "__")[[1]][2]
    1323. tax_table(ps.ng.tax_most_)[143,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[143,"Species"], "__")[[1]][2]
    1324. tax_table(ps.ng.tax_most_)[144,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[144,"Domain"], "__")[[1]][2]
    1325. tax_table(ps.ng.tax_most_)[144,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[144,"Phylum"], "__")[[1]][2]
    1326. tax_table(ps.ng.tax_most_)[144,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[144,"Class"], "__")[[1]][2]
    1327. tax_table(ps.ng.tax_most_)[144,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[144,"Order"], "__")[[1]][2]
    1328. tax_table(ps.ng.tax_most_)[144,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[144,"Family"], "__")[[1]][2]
    1329. tax_table(ps.ng.tax_most_)[144,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[144,"Genus"], "__")[[1]][2]
    1330. tax_table(ps.ng.tax_most_)[144,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[144,"Species"], "__")[[1]][2]
    1331. tax_table(ps.ng.tax_most_)[145,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[145,"Domain"], "__")[[1]][2]
    1332. tax_table(ps.ng.tax_most_)[145,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[145,"Phylum"], "__")[[1]][2]
    1333. tax_table(ps.ng.tax_most_)[145,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[145,"Class"], "__")[[1]][2]
    1334. tax_table(ps.ng.tax_most_)[145,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[145,"Order"], "__")[[1]][2]
    1335. tax_table(ps.ng.tax_most_)[145,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[145,"Family"], "__")[[1]][2]
    1336. tax_table(ps.ng.tax_most_)[145,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[145,"Genus"], "__")[[1]][2]
    1337. tax_table(ps.ng.tax_most_)[145,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[145,"Species"], "__")[[1]][2]
    1338. tax_table(ps.ng.tax_most_)[146,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[146,"Domain"], "__")[[1]][2]
    1339. tax_table(ps.ng.tax_most_)[146,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[146,"Phylum"], "__")[[1]][2]
    1340. tax_table(ps.ng.tax_most_)[146,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[146,"Class"], "__")[[1]][2]
    1341. tax_table(ps.ng.tax_most_)[146,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[146,"Order"], "__")[[1]][2]
    1342. tax_table(ps.ng.tax_most_)[146,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[146,"Family"], "__")[[1]][2]
    1343. tax_table(ps.ng.tax_most_)[146,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[146,"Genus"], "__")[[1]][2]
    1344. tax_table(ps.ng.tax_most_)[146,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[146,"Species"], "__")[[1]][2]
    1345. tax_table(ps.ng.tax_most_)[147,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[147,"Domain"], "__")[[1]][2]
    1346. tax_table(ps.ng.tax_most_)[147,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[147,"Phylum"], "__")[[1]][2]
    1347. tax_table(ps.ng.tax_most_)[147,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[147,"Class"], "__")[[1]][2]
    1348. tax_table(ps.ng.tax_most_)[147,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[147,"Order"], "__")[[1]][2]
    1349. tax_table(ps.ng.tax_most_)[147,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[147,"Family"], "__")[[1]][2]
    1350. tax_table(ps.ng.tax_most_)[147,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[147,"Genus"], "__")[[1]][2]
    1351. tax_table(ps.ng.tax_most_)[147,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[147,"Species"], "__")[[1]][2]
    1352. tax_table(ps.ng.tax_most_)[148,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[148,"Domain"], "__")[[1]][2]
    1353. tax_table(ps.ng.tax_most_)[148,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[148,"Phylum"], "__")[[1]][2]
    1354. tax_table(ps.ng.tax_most_)[148,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[148,"Class"], "__")[[1]][2]
    1355. tax_table(ps.ng.tax_most_)[148,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[148,"Order"], "__")[[1]][2]
    1356. tax_table(ps.ng.tax_most_)[148,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[148,"Family"], "__")[[1]][2]
    1357. tax_table(ps.ng.tax_most_)[148,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[148,"Genus"], "__")[[1]][2]
    1358. tax_table(ps.ng.tax_most_)[148,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[148,"Species"], "__")[[1]][2]
    1359. tax_table(ps.ng.tax_most_)[149,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[149,"Domain"], "__")[[1]][2]
    1360. tax_table(ps.ng.tax_most_)[149,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[149,"Phylum"], "__")[[1]][2]
    1361. tax_table(ps.ng.tax_most_)[149,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[149,"Class"], "__")[[1]][2]
    1362. tax_table(ps.ng.tax_most_)[149,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[149,"Order"], "__")[[1]][2]
    1363. tax_table(ps.ng.tax_most_)[149,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[149,"Family"], "__")[[1]][2]
    1364. tax_table(ps.ng.tax_most_)[149,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[149,"Genus"], "__")[[1]][2]
    1365. tax_table(ps.ng.tax_most_)[149,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[149,"Species"], "__")[[1]][2]
    1366. tax_table(ps.ng.tax_most_)[150,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[150,"Domain"], "__")[[1]][2]
    1367. tax_table(ps.ng.tax_most_)[150,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[150,"Phylum"], "__")[[1]][2]
    1368. tax_table(ps.ng.tax_most_)[150,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[150,"Class"], "__")[[1]][2]
    1369. tax_table(ps.ng.tax_most_)[150,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[150,"Order"], "__")[[1]][2]
    1370. tax_table(ps.ng.tax_most_)[150,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[150,"Family"], "__")[[1]][2]
    1371. tax_table(ps.ng.tax_most_)[150,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[150,"Genus"], "__")[[1]][2]
    1372. tax_table(ps.ng.tax_most_)[150,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[150,"Species"], "__")[[1]][2]
    1373. tax_table(ps.ng.tax_most_)[151,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[151,"Domain"], "__")[[1]][2]
    1374. tax_table(ps.ng.tax_most_)[151,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[151,"Phylum"], "__")[[1]][2]
    1375. tax_table(ps.ng.tax_most_)[151,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[151,"Class"], "__")[[1]][2]
    1376. tax_table(ps.ng.tax_most_)[151,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[151,"Order"], "__")[[1]][2]
    1377. tax_table(ps.ng.tax_most_)[151,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[151,"Family"], "__")[[1]][2]
    1378. tax_table(ps.ng.tax_most_)[151,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[151,"Genus"], "__")[[1]][2]
    1379. tax_table(ps.ng.tax_most_)[151,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[151,"Species"], "__")[[1]][2]
    1380. tax_table(ps.ng.tax_most_)[152,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[152,"Domain"], "__")[[1]][2]
    1381. tax_table(ps.ng.tax_most_)[152,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[152,"Phylum"], "__")[[1]][2]
    1382. tax_table(ps.ng.tax_most_)[152,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[152,"Class"], "__")[[1]][2]
    1383. tax_table(ps.ng.tax_most_)[152,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[152,"Order"], "__")[[1]][2]
    1384. tax_table(ps.ng.tax_most_)[152,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[152,"Family"], "__")[[1]][2]
    1385. tax_table(ps.ng.tax_most_)[152,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[152,"Genus"], "__")[[1]][2]
    1386. tax_table(ps.ng.tax_most_)[152,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[152,"Species"], "__")[[1]][2]
    1387. tax_table(ps.ng.tax_most_)[153,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[153,"Domain"], "__")[[1]][2]
    1388. tax_table(ps.ng.tax_most_)[153,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[153,"Phylum"], "__")[[1]][2]
    1389. tax_table(ps.ng.tax_most_)[153,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[153,"Class"], "__")[[1]][2]
    1390. tax_table(ps.ng.tax_most_)[153,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[153,"Order"], "__")[[1]][2]
    1391. tax_table(ps.ng.tax_most_)[153,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[153,"Family"], "__")[[1]][2]
    1392. tax_table(ps.ng.tax_most_)[153,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[153,"Genus"], "__")[[1]][2]
    1393. tax_table(ps.ng.tax_most_)[153,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[153,"Species"], "__")[[1]][2]
    1394. tax_table(ps.ng.tax_most_)[154,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[154,"Domain"], "__")[[1]][2]
    1395. tax_table(ps.ng.tax_most_)[154,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[154,"Phylum"], "__")[[1]][2]
    1396. tax_table(ps.ng.tax_most_)[154,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[154,"Class"], "__")[[1]][2]
    1397. tax_table(ps.ng.tax_most_)[154,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[154,"Order"], "__")[[1]][2]
    1398. tax_table(ps.ng.tax_most_)[154,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[154,"Family"], "__")[[1]][2]
    1399. tax_table(ps.ng.tax_most_)[154,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[154,"Genus"], "__")[[1]][2]
    1400. tax_table(ps.ng.tax_most_)[154,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[154,"Species"], "__")[[1]][2]
    1401. tax_table(ps.ng.tax_most_)[155,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[155,"Domain"], "__")[[1]][2]
    1402. tax_table(ps.ng.tax_most_)[155,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[155,"Phylum"], "__")[[1]][2]
    1403. tax_table(ps.ng.tax_most_)[155,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[155,"Class"], "__")[[1]][2]
    1404. tax_table(ps.ng.tax_most_)[155,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[155,"Order"], "__")[[1]][2]
    1405. tax_table(ps.ng.tax_most_)[155,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[155,"Family"], "__")[[1]][2]
    1406. tax_table(ps.ng.tax_most_)[155,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[155,"Genus"], "__")[[1]][2]
    1407. tax_table(ps.ng.tax_most_)[155,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[155,"Species"], "__")[[1]][2]
    1408. tax_table(ps.ng.tax_most_)[156,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[156,"Domain"], "__")[[1]][2]
    1409. tax_table(ps.ng.tax_most_)[156,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[156,"Phylum"], "__")[[1]][2]
    1410. tax_table(ps.ng.tax_most_)[156,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[156,"Class"], "__")[[1]][2]
    1411. tax_table(ps.ng.tax_most_)[156,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[156,"Order"], "__")[[1]][2]
    1412. tax_table(ps.ng.tax_most_)[156,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[156,"Family"], "__")[[1]][2]
    1413. tax_table(ps.ng.tax_most_)[156,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[156,"Genus"], "__")[[1]][2]
    1414. tax_table(ps.ng.tax_most_)[156,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[156,"Species"], "__")[[1]][2]
    1415. tax_table(ps.ng.tax_most_)[157,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[157,"Domain"], "__")[[1]][2]
    1416. tax_table(ps.ng.tax_most_)[157,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[157,"Phylum"], "__")[[1]][2]
    1417. tax_table(ps.ng.tax_most_)[157,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[157,"Class"], "__")[[1]][2]
    1418. tax_table(ps.ng.tax_most_)[157,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[157,"Order"], "__")[[1]][2]
    1419. tax_table(ps.ng.tax_most_)[157,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[157,"Family"], "__")[[1]][2]
    1420. tax_table(ps.ng.tax_most_)[157,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[157,"Genus"], "__")[[1]][2]
    1421. tax_table(ps.ng.tax_most_)[157,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[157,"Species"], "__")[[1]][2]
    1422. tax_table(ps.ng.tax_most_)[158,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[158,"Domain"], "__")[[1]][2]
    1423. tax_table(ps.ng.tax_most_)[158,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[158,"Phylum"], "__")[[1]][2]
    1424. tax_table(ps.ng.tax_most_)[158,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[158,"Class"], "__")[[1]][2]
    1425. tax_table(ps.ng.tax_most_)[158,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[158,"Order"], "__")[[1]][2]
    1426. tax_table(ps.ng.tax_most_)[158,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[158,"Family"], "__")[[1]][2]
    1427. tax_table(ps.ng.tax_most_)[158,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[158,"Genus"], "__")[[1]][2]
    1428. tax_table(ps.ng.tax_most_)[158,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[158,"Species"], "__")[[1]][2]
    1429. tax_table(ps.ng.tax_most_)[159,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[159,"Domain"], "__")[[1]][2]
    1430. tax_table(ps.ng.tax_most_)[159,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[159,"Phylum"], "__")[[1]][2]
    1431. tax_table(ps.ng.tax_most_)[159,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[159,"Class"], "__")[[1]][2]
    1432. tax_table(ps.ng.tax_most_)[159,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[159,"Order"], "__")[[1]][2]
    1433. tax_table(ps.ng.tax_most_)[159,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[159,"Family"], "__")[[1]][2]
    1434. tax_table(ps.ng.tax_most_)[159,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[159,"Genus"], "__")[[1]][2]
    1435. tax_table(ps.ng.tax_most_)[159,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[159,"Species"], "__")[[1]][2]
    1436. tax_table(ps.ng.tax_most_)[160,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[160,"Domain"], "__")[[1]][2]
    1437. tax_table(ps.ng.tax_most_)[160,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[160,"Phylum"], "__")[[1]][2]
    1438. tax_table(ps.ng.tax_most_)[160,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[160,"Class"], "__")[[1]][2]
    1439. tax_table(ps.ng.tax_most_)[160,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[160,"Order"], "__")[[1]][2]
    1440. tax_table(ps.ng.tax_most_)[160,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[160,"Family"], "__")[[1]][2]
    1441. tax_table(ps.ng.tax_most_)[160,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[160,"Genus"], "__")[[1]][2]
    1442. tax_table(ps.ng.tax_most_)[160,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[160,"Species"], "__")[[1]][2]
    1443. tax_table(ps.ng.tax_most_)[161,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[161,"Domain"], "__")[[1]][2]
    1444. tax_table(ps.ng.tax_most_)[161,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[161,"Phylum"], "__")[[1]][2]
    1445. tax_table(ps.ng.tax_most_)[161,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[161,"Class"], "__")[[1]][2]
    1446. tax_table(ps.ng.tax_most_)[161,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[161,"Order"], "__")[[1]][2]
    1447. tax_table(ps.ng.tax_most_)[161,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[161,"Family"], "__")[[1]][2]
    1448. tax_table(ps.ng.tax_most_)[161,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[161,"Genus"], "__")[[1]][2]
    1449. tax_table(ps.ng.tax_most_)[161,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[161,"Species"], "__")[[1]][2]
    1450. tax_table(ps.ng.tax_most_)[162,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[162,"Domain"], "__")[[1]][2]
    1451. tax_table(ps.ng.tax_most_)[162,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[162,"Phylum"], "__")[[1]][2]
    1452. tax_table(ps.ng.tax_most_)[162,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[162,"Class"], "__")[[1]][2]
    1453. tax_table(ps.ng.tax_most_)[162,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[162,"Order"], "__")[[1]][2]
    1454. tax_table(ps.ng.tax_most_)[162,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[162,"Family"], "__")[[1]][2]
    1455. tax_table(ps.ng.tax_most_)[162,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[162,"Genus"], "__")[[1]][2]
    1456. tax_table(ps.ng.tax_most_)[162,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[162,"Species"], "__")[[1]][2]
    1457. tax_table(ps.ng.tax_most_)[163,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[163,"Domain"], "__")[[1]][2]
    1458. tax_table(ps.ng.tax_most_)[163,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[163,"Phylum"], "__")[[1]][2]
    1459. tax_table(ps.ng.tax_most_)[163,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[163,"Class"], "__")[[1]][2]
    1460. tax_table(ps.ng.tax_most_)[163,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[163,"Order"], "__")[[1]][2]
    1461. tax_table(ps.ng.tax_most_)[163,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[163,"Family"], "__")[[1]][2]
    1462. tax_table(ps.ng.tax_most_)[163,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[163,"Genus"], "__")[[1]][2]
    1463. tax_table(ps.ng.tax_most_)[163,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[163,"Species"], "__")[[1]][2]
    1464. tax_table(ps.ng.tax_most_)[164,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[164,"Domain"], "__")[[1]][2]
    1465. tax_table(ps.ng.tax_most_)[164,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[164,"Phylum"], "__")[[1]][2]
    1466. tax_table(ps.ng.tax_most_)[164,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[164,"Class"], "__")[[1]][2]
    1467. tax_table(ps.ng.tax_most_)[164,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[164,"Order"], "__")[[1]][2]
    1468. tax_table(ps.ng.tax_most_)[164,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[164,"Family"], "__")[[1]][2]
    1469. tax_table(ps.ng.tax_most_)[164,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[164,"Genus"], "__")[[1]][2]
    1470. tax_table(ps.ng.tax_most_)[164,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[164,"Species"], "__")[[1]][2]
    1471. tax_table(ps.ng.tax_most_)[165,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[165,"Domain"], "__")[[1]][2]
    1472. tax_table(ps.ng.tax_most_)[165,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[165,"Phylum"], "__")[[1]][2]
    1473. tax_table(ps.ng.tax_most_)[165,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[165,"Class"], "__")[[1]][2]
    1474. tax_table(ps.ng.tax_most_)[165,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[165,"Order"], "__")[[1]][2]
    1475. tax_table(ps.ng.tax_most_)[165,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[165,"Family"], "__")[[1]][2]
    1476. tax_table(ps.ng.tax_most_)[165,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[165,"Genus"], "__")[[1]][2]
    1477. tax_table(ps.ng.tax_most_)[165,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[165,"Species"], "__")[[1]][2]
    1478. tax_table(ps.ng.tax_most_)[166,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[166,"Domain"], "__")[[1]][2]
    1479. tax_table(ps.ng.tax_most_)[166,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[166,"Phylum"], "__")[[1]][2]
    1480. tax_table(ps.ng.tax_most_)[166,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[166,"Class"], "__")[[1]][2]
    1481. tax_table(ps.ng.tax_most_)[166,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[166,"Order"], "__")[[1]][2]
    1482. tax_table(ps.ng.tax_most_)[166,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[166,"Family"], "__")[[1]][2]
    1483. tax_table(ps.ng.tax_most_)[166,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[166,"Genus"], "__")[[1]][2]
    1484. tax_table(ps.ng.tax_most_)[166,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[166,"Species"], "__")[[1]][2]
    1485. tax_table(ps.ng.tax_most_)[167,"Domain"] <- str_split(tax_table(ps.ng.tax_most_)[167,"Domain"], "__")[[1]][2]
    1486. tax_table(ps.ng.tax_most_)[167,"Phylum"] <- str_split(tax_table(ps.ng.tax_most_)[167,"Phylum"], "__")[[1]][2]
    1487. tax_table(ps.ng.tax_most_)[167,"Class"] <- str_split(tax_table(ps.ng.tax_most_)[167,"Class"], "__")[[1]][2]
    1488. tax_table(ps.ng.tax_most_)[167,"Order"] <- str_split(tax_table(ps.ng.tax_most_)[167,"Order"], "__")[[1]][2]
    1489. tax_table(ps.ng.tax_most_)[167,"Family"] <- str_split(tax_table(ps.ng.tax_most_)[167,"Family"], "__")[[1]][2]
    1490. tax_table(ps.ng.tax_most_)[167,"Genus"] <- str_split(tax_table(ps.ng.tax_most_)[167,"Genus"], "__")[[1]][2]
    1491. tax_table(ps.ng.tax_most_)[167,"Species"] <- str_split(tax_table(ps.ng.tax_most_)[167,"Species"], "__")[[1]][2]
    1492. ```
    1493. ```{r, echo=TRUE, warning=FALSE}
    1494. #aes(color="Phylum", fill="Phylum") --> aes()
    1495. #ggplot(data=data, aes(x=Sample, y=Abundance, fill=Phylum))
    1496. plot_bar(ps.ng.tax_most_, fill="Phylum") + geom_bar(aes(), stat="identity", position="stack") +
    1497. scale_fill_manual(values = c("darkblue", "darkgoldenrod1", "darkseagreen", "darkorchid", "darkolivegreen1", "lightskyblue", "darkgreen", "deeppink", "khaki2", "firebrick", "brown1", "darkorange1", "cyan1", "royalblue4", "darksalmon", "darkblue","royalblue4", "dodgerblue3", "steelblue1", "lightskyblue", "darkseagreen", "darkgoldenrod1", "darkseagreen", "darkorchid", "darkolivegreen1", "brown1", "darkorange1", "cyan1", "darkgrey")) + theme(axis.text = element_text(size = 5, colour="black")) + theme(legend.position="bottom") + guides(fill=guide_legend(nrow=2)) #6 instead of theme.size
    1498. ```
    1499. ```{r, echo=FALSE, warning=FALSE}
    1500. #png("abc.png")
    1501. #knitr::include_graphics("./Phyloseq_files/figure-html/unnamed-chunk-7-1.png")
    1502. #dev.off()
    1503. ```
    1504. \pagebreak
    1505. Regroup together pre vs post stroke samples and normalize number of reads in each group using median sequencing depth.
    1506. ```{r, echo=TRUE, warning=FALSE}
    1507. ps.ng.tax_most_pre_post_stroke <- merge_samples(ps.ng.tax_most_, "pre_post_stroke")
    1508. ps.ng.tax_most_pre_post_stroke_ = transform_sample_counts(ps.ng.tax_most_pre_post_stroke, function(x) x / sum(x))
    1509. #plot_bar(ps.ng.tax_most_SampleType_, fill = "Phylum") + geom_bar(aes(color=Phylum, fill=Phylum), stat="identity", position="stack")
    1510. plot_bar(ps.ng.tax_most_pre_post_stroke_, fill="Phylum") + geom_bar(aes(), stat="identity", position="stack") +
    1511. scale_fill_manual(values = c("darkblue", "darkgoldenrod1", "darkseagreen", "darkorchid", "darkolivegreen1", "lightskyblue", "darkgreen", "deeppink", "khaki2", "firebrick", "brown1", "darkorange1", "cyan1", "royalblue4", "darksalmon", "darkblue","royalblue4", "dodgerblue3", "steelblue1", "lightskyblue", "darkseagreen", "darkgoldenrod1", "darkseagreen", "darkorchid", "darkolivegreen1", "brown1", "darkorange1", "cyan1", "darkgrey")) + theme(axis.text = element_text(size = theme.size, colour="black"))
    1512. ```
    1513. \pagebreak
    1514. Use color according to phylum. Do separate panels Stroke and Sex_age.
    1515. ```{r, echo=TRUE, warning=FALSE}
    1516. ps.ng.tax_most_copied <- data.table::copy(ps.ng.tax_most_)
    1517. #FITTING6: regulate the bar height if it has replicates: 5+6+6+8+4+8+6+7=25+25=50
    1518. otu_table(ps.ng.tax_most_)[,c("1")] <- otu_table(ps.ng.tax_most_)[,c("1")]/5
    1519. otu_table(ps.ng.tax_most_)[,c("2")] <- otu_table(ps.ng.tax_most_)[,c("2")]/5
    1520. otu_table(ps.ng.tax_most_)[,c("5")] <- otu_table(ps.ng.tax_most_)[,c("5")]/5
    1521. otu_table(ps.ng.tax_most_)[,c("6")] <- otu_table(ps.ng.tax_most_)[,c("6")]/5
    1522. otu_table(ps.ng.tax_most_)[,c("7")] <- otu_table(ps.ng.tax_most_)[,c("7")]/5
    1523. otu_table(ps.ng.tax_most_)[,c("8")] <- otu_table(ps.ng.tax_most_)[,c("8")]/6
    1524. otu_table(ps.ng.tax_most_)[,c("9")] <- otu_table(ps.ng.tax_most_)[,c("9")]/6
    1525. otu_table(ps.ng.tax_most_)[,c("10")] <- otu_table(ps.ng.tax_most_)[,c("10")]/6
    1526. otu_table(ps.ng.tax_most_)[,c("12")] <- otu_table(ps.ng.tax_most_)[,c("12")]/6
    1527. otu_table(ps.ng.tax_most_)[,c("13")] <- otu_table(ps.ng.tax_most_)[,c("13")]/6
    1528. otu_table(ps.ng.tax_most_)[,c("14")] <- otu_table(ps.ng.tax_most_)[,c("14")]/6
    1529. otu_table(ps.ng.tax_most_)[,c("15")] <- otu_table(ps.ng.tax_most_)[,c("15")]/6
    1530. otu_table(ps.ng.tax_most_)[,c("16")] <- otu_table(ps.ng.tax_most_)[,c("16")]/6
    1531. otu_table(ps.ng.tax_most_)[,c("17")] <- otu_table(ps.ng.tax_most_)[,c("17")]/6
    1532. otu_table(ps.ng.tax_most_)[,c("18")] <- otu_table(ps.ng.tax_most_)[,c("18")]/6
    1533. otu_table(ps.ng.tax_most_)[,c("19")] <- otu_table(ps.ng.tax_most_)[,c("19")]/6
    1534. otu_table(ps.ng.tax_most_)[,c("20")] <- otu_table(ps.ng.tax_most_)[,c("20")]/6
    1535. otu_table(ps.ng.tax_most_)[,c("21")] <- otu_table(ps.ng.tax_most_)[,c("21")]/8
    1536. otu_table(ps.ng.tax_most_)[,c("22")] <- otu_table(ps.ng.tax_most_)[,c("22")]/8
    1537. otu_table(ps.ng.tax_most_)[,c("23")] <- otu_table(ps.ng.tax_most_)[,c("23")]/8
    1538. otu_table(ps.ng.tax_most_)[,c("24")] <- otu_table(ps.ng.tax_most_)[,c("24")]/8
    1539. otu_table(ps.ng.tax_most_)[,c("25")] <- otu_table(ps.ng.tax_most_)[,c("25")]/8
    1540. otu_table(ps.ng.tax_most_)[,c("26")] <- otu_table(ps.ng.tax_most_)[,c("26")]/8
    1541. otu_table(ps.ng.tax_most_)[,c("27")] <- otu_table(ps.ng.tax_most_)[,c("27")]/8
    1542. otu_table(ps.ng.tax_most_)[,c("28")] <- otu_table(ps.ng.tax_most_)[,c("28")]/8
    1543. otu_table(ps.ng.tax_most_)[,c("29")] <- otu_table(ps.ng.tax_most_)[,c("29")]/4
    1544. otu_table(ps.ng.tax_most_)[,c("30")] <- otu_table(ps.ng.tax_most_)[,c("30")]/4
    1545. otu_table(ps.ng.tax_most_)[,c("31")] <- otu_table(ps.ng.tax_most_)[,c("31")]/4
    1546. otu_table(ps.ng.tax_most_)[,c("32")] <- otu_table(ps.ng.tax_most_)[,c("32")]/4
    1547. otu_table(ps.ng.tax_most_)[,c("33")] <- otu_table(ps.ng.tax_most_)[,c("33")]/8
    1548. otu_table(ps.ng.tax_most_)[,c("34")] <- otu_table(ps.ng.tax_most_)[,c("34")]/8
    1549. otu_table(ps.ng.tax_most_)[,c("35")] <- otu_table(ps.ng.tax_most_)[,c("35")]/8
    1550. otu_table(ps.ng.tax_most_)[,c("36")] <- otu_table(ps.ng.tax_most_)[,c("36")]/8
    1551. otu_table(ps.ng.tax_most_)[,c("37")] <- otu_table(ps.ng.tax_most_)[,c("37")]/8
    1552. otu_table(ps.ng.tax_most_)[,c("38")] <- otu_table(ps.ng.tax_most_)[,c("38")]/8
    1553. otu_table(ps.ng.tax_most_)[,c("39")] <- otu_table(ps.ng.tax_most_)[,c("39")]/8
    1554. otu_table(ps.ng.tax_most_)[,c("51")] <- otu_table(ps.ng.tax_most_)[,c("51")]/8
    1555. otu_table(ps.ng.tax_most_)[,c("40")] <- otu_table(ps.ng.tax_most_)[,c("40")]/6
    1556. otu_table(ps.ng.tax_most_)[,c("41")] <- otu_table(ps.ng.tax_most_)[,c("41")]/6
    1557. otu_table(ps.ng.tax_most_)[,c("42")] <- otu_table(ps.ng.tax_most_)[,c("42")]/6
    1558. otu_table(ps.ng.tax_most_)[,c("43")] <- otu_table(ps.ng.tax_most_)[,c("43")]/6
    1559. otu_table(ps.ng.tax_most_)[,c("44")] <- otu_table(ps.ng.tax_most_)[,c("44")]/6
    1560. otu_table(ps.ng.tax_most_)[,c("46")] <- otu_table(ps.ng.tax_most_)[,c("46")]/6
    1561. otu_table(ps.ng.tax_most_)[,c("47")] <- otu_table(ps.ng.tax_most_)[,c("47")]/7
    1562. otu_table(ps.ng.tax_most_)[,c("48")] <- otu_table(ps.ng.tax_most_)[,c("48")]/7
    1563. otu_table(ps.ng.tax_most_)[,c("49")] <- otu_table(ps.ng.tax_most_)[,c("49")]/7
    1564. otu_table(ps.ng.tax_most_)[,c("50")] <- otu_table(ps.ng.tax_most_)[,c("50")]/7
    1565. otu_table(ps.ng.tax_most_)[,c("52")] <- otu_table(ps.ng.tax_most_)[,c("52")]/7
    1566. otu_table(ps.ng.tax_most_)[,c("53")] <- otu_table(ps.ng.tax_most_)[,c("53")]/7
    1567. otu_table(ps.ng.tax_most_)[,c("55")] <- otu_table(ps.ng.tax_most_)[,c("55")]/7
    1568. #plot_bar(ps.ng.tax_most_swab_, x="Phylum", fill = "Phylum", facet_grid = Patient~RoundDay) + geom_bar(aes(color=Phylum, fill=Phylum), stat="identity", position="stack") + theme(axis.text = element_text(size = theme.size, colour="black"))
    1569. plot_bar(ps.ng.tax_most_, x="Phylum", fill="Phylum", facet_grid = pre_post_stroke~Sex_age) + geom_bar(aes(), stat="identity", position="stack") +
    1570. scale_fill_manual(values = c("darkblue", "darkgoldenrod1", "darkseagreen", "darkorchid", "darkolivegreen1", "lightskyblue", "darkgreen", "deeppink", "khaki2", "firebrick", "brown1", "darkorange1", "cyan1", "royalblue4", "darksalmon", "darkblue","royalblue4", "dodgerblue3", "steelblue1", "lightskyblue", "darkseagreen", "darkgoldenrod1", "darkseagreen", "darkorchid", "darkolivegreen1", "brown1", "darkorange1", "cyan1", "darkgrey")) + theme(axis.text = element_text(size = 5, colour="black"), axis.text.x=element_blank(), axis.ticks=element_blank()) + theme(legend.position="bottom") + guides(fill=guide_legend(nrow=2))
    1571. ```
    1572. ```{r, echo=FALSE, warning=FALSE}
    1573. #knitr::include_graphics("./Phyloseq_files/figure-html/unnamed-chunk-10-1.png")
    1574. #> tax_table(carbom)
    1575. #Taxonomy Table: [205 taxa by 7 taxonomic ranks]:
    1576. # Domain Supergroup Division Class
    1577. #Otu001 "Eukaryota" "Archaeplastida" "Chlorophyta" "Mamiellophyceae"
    1578. # Order Family Genus
    1579. #Otu001 "Mamiellales" "Bathycoccaceae" "Ostreococcus"
    1580. #sample_data(ps.ng.tax)
    1581. ```
    1582. ## Bar plots in class level
    1583. ```{r, echo=TRUE, warning=FALSE}
    1584. plot_bar(ps.ng.tax_most_copied, fill="Class") + geom_bar(aes(), stat="identity", position="stack") +
    1585. scale_fill_manual(values = c("darkblue", "darkgoldenrod1", "darkseagreen", "darkorchid", "darkolivegreen1", "lightskyblue", "darkgreen", "deeppink", "khaki2", "firebrick", "brown1", "darkorange1", "cyan1", "royalblue4", "darksalmon", "darkblue","royalblue4", "dodgerblue3", "steelblue1", "lightskyblue", "darkseagreen", "darkgoldenrod1", "darkseagreen", "darkorchid", "darkolivegreen1", "brown1", "darkorange1", "cyan1", "darkgrey")) + theme(axis.text = element_text(size = 5, colour="black")) + theme(legend.position="bottom") + guides(fill=guide_legend(nrow=3))
    1586. ```
    1587. Regroup together pre vs post stroke samples and normalize number of reads in each group using median sequencing depth.
    1588. ```{r, echo=TRUE, warning=FALSE}
    1589. plot_bar(ps.ng.tax_most_pre_post_stroke_, fill="Class") + geom_bar(aes(), stat="identity", position="stack") +
    1590. scale_fill_manual(values = c("darkblue", "darkgoldenrod1", "darkseagreen", "darkorchid", "darkolivegreen1", "lightskyblue", "darkgreen", "deeppink", "khaki2", "firebrick", "brown1", "darkorange1", "cyan1", "royalblue4", "darksalmon", "darkblue","royalblue4", "dodgerblue3", "steelblue1", "lightskyblue", "darkseagreen", "darkgoldenrod1", "darkseagreen", "darkorchid", "darkolivegreen1", "brown1", "darkorange1", "cyan1", "darkgrey")) + theme(axis.text = element_text(size = theme.size, colour="black"))
    1591. ```
    1592. \pagebreak
    1593. Use color according to class. Do separate panels Stroke and Sex_age.
    1594. ```{r, echo=TRUE, warning=FALSE}
    1595. #-- If existing replicates, to be processed as follows --
    1596. plot_bar(ps.ng.tax_most_, x="Class", fill="Class", facet_grid = pre_post_stroke~Sex_age) + geom_bar(aes(), stat="identity", position="stack") +
    1597. scale_fill_manual(values = c("darkblue", "darkgoldenrod1", "darkseagreen", "darkorchid", "darkolivegreen1", "lightskyblue", "darkgreen", "deeppink", "khaki2", "firebrick", "brown1", "darkorange1", "cyan1", "royalblue4", "darksalmon", "darkblue","royalblue4", "dodgerblue3", "steelblue1", "lightskyblue", "darkseagreen", "darkgoldenrod1", "darkseagreen", "darkorchid", "darkolivegreen1", "brown1", "darkorange1", "cyan1", "darkgrey")) + theme(axis.text = element_text(size = 5, colour="black"), axis.text.x=element_blank(), axis.ticks=element_blank()) + theme(legend.position="bottom") + guides(fill=guide_legend(nrow=3))
    1598. ```
    1599. ## Bar plots in order level
    1600. ```{r, echo=TRUE, warning=FALSE}
    1601. plot_bar(ps.ng.tax_most_copied, fill="Order") + geom_bar(aes(), stat="identity", position="stack") +
    1602. scale_fill_manual(values = c("darkblue", "darkgoldenrod1", "darkseagreen", "darkorchid", "darkolivegreen1", "lightskyblue", "darkgreen", "deeppink", "khaki2", "firebrick", "brown1", "darkorange1", "cyan1", "royalblue4", "darksalmon", "darkblue","royalblue4", "dodgerblue3", "steelblue1", "lightskyblue", "darkseagreen", "darkgoldenrod1", "darkseagreen", "darkorchid", "darkolivegreen1", "brown1", "darkorange1", "cyan1", "darkgrey")) + theme(axis.text = element_text(size = 5, colour="black")) + theme(legend.position="bottom") + guides(fill=guide_legend(nrow=4))
    1603. ```
    1604. Regroup together pre vs post stroke and normalize number of reads in each group using median sequencing depth.
    1605. ```{r, echo=TRUE, warning=FALSE}
    1606. plot_bar(ps.ng.tax_most_pre_post_stroke_, fill="Order") + geom_bar(aes(), stat="identity", position="stack") +
    1607. scale_fill_manual(values = c("darkblue", "darkgoldenrod1", "darkseagreen", "darkorchid", "darkolivegreen1", "lightskyblue", "darkgreen", "deeppink", "khaki2", "firebrick", "brown1", "darkorange1", "cyan1", "royalblue4", "darksalmon", "darkblue","royalblue4", "dodgerblue3", "steelblue1", "lightskyblue", "darkseagreen", "darkgoldenrod1", "darkseagreen", "darkorchid", "darkolivegreen1", "brown1", "darkorange1", "cyan1", "darkgrey")) + theme(axis.text = element_text(size = theme.size, colour="black")) + theme(legend.position="bottom") + guides(fill=guide_legend(nrow=4))
    1608. ```
    1609. \pagebreak
    1610. Use color according to order. Do separate panels Stroke and Sex_age.
    1611. ```{r, echo=TRUE, warning=FALSE}
    1612. #FITTING7: regulate the bar height if it has replicates
    1613. plot_bar(ps.ng.tax_most_, x="Order", fill="Order", facet_grid = pre_post_stroke~Sex_age) + geom_bar(aes(), stat="identity", position="stack") +
    1614. scale_fill_manual(values = c("darkblue", "darkgoldenrod1", "darkseagreen", "darkorchid", "darkolivegreen1", "lightskyblue", "darkgreen", "deeppink", "khaki2", "firebrick", "brown1", "darkorange1", "cyan1", "royalblue4", "darksalmon", "darkblue","royalblue4", "dodgerblue3", "steelblue1", "lightskyblue", "darkseagreen", "darkgoldenrod1", "darkseagreen", "darkorchid", "darkolivegreen1", "brown1", "darkorange1", "cyan1", "darkgrey")) + theme(axis.text = element_text(size = 5, colour="black"), axis.text.x=element_blank(), axis.ticks=element_blank()) + theme(legend.position="bottom") + guides(fill=guide_legend(nrow=4))
    1615. ```
    1616. ## Bar plots in family level
    1617. ```{r, echo=TRUE, warning=FALSE}
    1618. plot_bar(ps.ng.tax_most_copied, fill="Family") + geom_bar(aes(), stat="identity", position="stack") +
    1619. scale_fill_manual(values = c("#FF0000", "#000000", "#0000FF", "#C0C0C0", "#FFFFFF", "#FFFF00", "#00FFFF", "#FFA500", "#00FF00", "#808080", "#FF00FF", "#800080", "#FDD017", "#0000A0", "#3BB9FF", "#008000", "#800000", "#ADD8E6", "#F778A1", "#800517", "#736F6E", "#F52887", "#C11B17", "#5CB3FF", "#A52A2A", "#FF8040", "#2B60DE", "#736AFF", "#1589FF", "#98AFC7", "#8D38C9", "#307D7E", "#F6358A", "#151B54", "#6D7B8D", "#FDEEF4", "#FF0080", "#F88017", "#2554C7", "#FFF8C6", "#D4A017", "#306EFF", "#151B8D", "#9E7BFF", "#EAC117", "#E0FFFF", "#15317E", "#6C2DC7", "#FBB917", "#FCDFFF", "#15317E", "#254117", "#FAAFBE", "#357EC7")) + theme(axis.text = element_text(size = 5, colour="black")) + theme(legend.position="bottom") + guides(fill=guide_legend(nrow=8))
    1620. ```
    1621. Regroup together pre vs post stroke samples and normalize number of reads in each group using median sequencing depth.
    1622. ```{r, echo=TRUE, warning=FALSE}
    1623. plot_bar(ps.ng.tax_most_pre_post_stroke_, fill="Family") + geom_bar(aes(), stat="identity", position="stack") +
    1624. scale_fill_manual(values = c("#FF0000", "#000000", "#0000FF", "#C0C0C0", "#FFFFFF", "#FFFF00", "#00FFFF", "#FFA500", "#00FF00", "#808080", "#FF00FF", "#800080", "#FDD017", "#0000A0", "#3BB9FF", "#008000", "#800000", "#ADD8E6", "#F778A1", "#800517", "#736F6E", "#F52887", "#C11B17", "#5CB3FF", "#A52A2A", "#FF8040", "#2B60DE", "#736AFF", "#1589FF", "#98AFC7", "#8D38C9", "#307D7E", "#F6358A", "#151B54", "#6D7B8D", "#FDEEF4", "#FF0080", "#F88017", "#2554C7", "#FFF8C6", "#D4A017", "#306EFF", "#151B8D", "#9E7BFF", "#EAC117", "#E0FFFF", "#15317E", "#6C2DC7", "#FBB917", "#FCDFFF", "#15317E", "#254117", "#FAAFBE", "#357EC7")) + theme(axis.text = element_text(size = theme.size, colour="black")) + theme(legend.position="bottom") + guides(fill=guide_legend(nrow=8))
    1625. ```
    1626. \pagebreak
    1627. Use color according to family. Do separate panels Stroke and Sex_age.
    1628. ```{r, echo=TRUE, warning=FALSE}
    1629. #-- If existing replicates, to be processed as follows --
    1630. plot_bar(ps.ng.tax_most_, x="Family", fill="Family", facet_grid = pre_post_stroke~Sex_age) + geom_bar(aes(), stat="identity", position="stack") +
    1631. scale_fill_manual(values = c("#FF0000", "#000000", "#0000FF", "#C0C0C0", "#FFFFFF", "#FFFF00", "#00FFFF", "#FFA500", "#00FF00", "#808080", "#FF00FF", "#800080", "#FDD017", "#0000A0", "#3BB9FF", "#008000", "#800000", "#ADD8E6", "#F778A1", "#800517", "#736F6E", "#F52887", "#C11B17", "#5CB3FF", "#A52A2A", "#FF8040", "#2B60DE", "#736AFF", "#1589FF", "#98AFC7", "#8D38C9", "#307D7E", "#F6358A", "#151B54", "#6D7B8D", "#FDEEF4", "#FF0080", "#F88017", "#2554C7", "#FFF8C6", "#D4A017", "#306EFF", "#151B8D", "#9E7BFF", "#EAC117", "#E0FFFF", "#15317E", "#6C2DC7", "#FBB917", "#FCDFFF", "#15317E", "#254117", "#FAAFBE", "#357EC7")) + theme(axis.text = element_text(size = 5, colour="black"), axis.text.x=element_blank(), axis.ticks=element_blank()) + theme(legend.position="bottom") + guides(fill=guide_legend(nrow=8))
    1632. ```
    1633. ```{r, echo=FALSE, warning=FALSE}
    1634. #MOVE_FROM_ABOVE: ## Bar plots in genus level
    1635. #MOVE_FROM_ABOVE: Regroup together pre vs post stroke samples and normalize number of reads in each group using median sequencing depth.
    1636. #plot_bar(ps.ng.tax_most_pre_post_stroke_, fill="Genus") + geom_bar(aes(), stat="identity", position="stack") +
    1637. #scale_fill_manual(values = c("#FF0000", "#000000", "#0000FF", "#C0C0C0", "#FFFFFF", "#FFFF00", "#00FFFF", "#FFA500", "#00FF00", "#808080", "#FF00FF", "#800080", "#FDD017", "#0000A0", "#3BB9FF", "#008000", "#800000", "#ADD8E6", "#F778A1", "#800517", "#736F6E", "#F52887", "#C11B17", "#5CB3FF", "#A52A2A", "#FF8040", "#2B60DE", "#736AFF", "#1589FF", "#98AFC7", "#8D38C9", "#307D7E", "#F6358A", "#151B54", "#6D7B8D", "#FDEEF4", "#FF0080", "#F88017", "#2554C7", "#FFF8C6", "#D4A017", "#306EFF", "#151B8D", "#9E7BFF", "#EAC117", "#E0FFFF", "#15317E", "#6C2DC7", "#FBB917", "#FCDFFF", "#15317E", "#254117", "#FAAFBE", "#357EC7")) + theme(axis.text = element_text(size = theme.size, colour="black")) + theme(legend.position="bottom")
    1638. ```
    1639. \pagebreak
    1640. ```{r, echo=FALSE, warning=FALSE}
    1641. #MOVE_FROM_ABOVE: Use color according to genus. Do separate panels Stroke and Sex_age.
    1642. ##-- If existing replicates, to be processed as follows --
    1643. #plot_bar(ps.ng.tax_most_, x="Genus", fill="Genus", facet_grid = pre_post_stroke~Sex_age) + geom_bar(aes(), stat="identity", position="stack") +
    1644. #scale_fill_manual(values = c("#FF0000", "#000000", "#0000FF", "#C0C0C0", "#FFFFFF", "#FFFF00", "#00FFFF", "#FFA500", "#00FF00", "#808080", "#FF00FF", "#800080", "#FDD017", "#0000A0", "#3BB9FF", "#008000", "#800000", "#ADD8E6", "#F778A1", "#800517", "#736F6E", "#F52887", "#C11B17", "#5CB3FF", "#A52A2A", "#FF8040", "#2B60DE", "#736AFF", "#1589FF", "#98AFC7", "#8D38C9", "#307D7E", "#F6358A", "#151B54", "#6D7B8D", "#FDEEF4", "#FF0080", "#F88017", "#2554C7", "#FFF8C6", "#D4A017", "#306EFF", "#151B8D", "#9E7BFF", "#EAC117", "#E0FFFF", "#15317E", "#6C2DC7", "#FBB917", "#FCDFFF", "#15317E", "#254117", "#FAAFBE", "#357EC7")) + theme(axis.text = element_text(size = 6, colour="black"), axis.text.x=element_blank(), axis.ticks=element_blank()) + theme(legend.position="bottom") + guides(fill=guide_legend(nrow=18))
    1645. ```
    1646. \pagebreak
    1647. # Alpha diversity
    1648. Plot Chao1 richness estimator, Observed OTUs, Shannon index, and Phylogenetic diversity.
    1649. Regroup together samples from the same group.
    1650. ```{r, echo=FALSE, warning=FALSE}
    1651. # using rarefied data
    1652. #FITTING2: CONSOLE:
    1653. #gunzip table_even4753.biom.gz
    1654. #alpha_diversity.py -i table_even42369.biom --metrics chao1,observed_otus,shannon,PD_whole_tree -o adiv_even.txt -t ../clustering/rep_set.tre
    1655. #gunzip table_even4753.biom.gz
    1656. #alpha_diversity.py -i table_even4753.biom --metrics chao1,observed_otus,shannon,PD_whole_tree -o adiv_even.txt -t ../clustering_stool/rep_set.tre
    1657. #gunzip table_even4753.biom.gz
    1658. #alpha_diversity.py -i table_even4753.biom --metrics chao1,observed_otus,shannon,PD_whole_tree -o adiv_even.txt -t ../clustering_swab/rep_set.tre
    1659. ```
    1660. ```{r, echo=TRUE, warning=FALSE}
    1661. hmp.div_qiime <- read.csv("adiv_even.txt", sep="\t")
    1662. colnames(hmp.div_qiime) <- c("sam_name", "chao1", "observed_otus", "shannon", "PD_whole_tree")
    1663. row.names(hmp.div_qiime) <- hmp.div_qiime$sam_name
    1664. div.df <- merge(hmp.div_qiime, hmp.meta, by = "sam_name")
    1665. div.df2 <- div.df[, c("Group", "chao1", "shannon", "observed_otus", "PD_whole_tree")]
    1666. colnames(div.df2) <- c("Group", "Chao-1", "Shannon", "OTU", "Phylogenetic Diversity")
    1667. #colnames(div.df2)
    1668. options(max.print=999999)
    1669. #27 H47 830.5000 5.008482 319 10.60177
    1670. #FITTING4: if occuring "Computation failed in `stat_signif()`:not enough 'y' observations"
    1671. #means: the patient H47 contains only one sample, it should be removed for the statistical p-values calculations.
    1672. #delete H47(1)
    1673. #div.df2 <- div.df2[-c(3), ]
    1674. #div.df2 <- div.df2[-c(55,54, 45,40,39,27,26,25,1), ]
    1675. knitr::kable(div.df2) %>% kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"))
    1676. #https://uc-r.github.io/t_test
    1677. #We can perform the test with t.test and transform our data and we can also perform the nonparametric test with the wilcox.test function.
    1678. stat.test.Shannon <- compare_means(
    1679. Shannon ~ Group, data = div.df2,
    1680. method = "t.test"
    1681. )
    1682. knitr::kable(stat.test.Shannon) %>% kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"))
    1683. div_df_melt <- reshape2::melt(div.df2)
    1684. #head(div_df_melt)
    1685. #https://plot.ly/r/box-plots/#horizontal-boxplot
    1686. #http://www.sthda.com/english/wiki/print.php?id=177
    1687. #https://rpkgs.datanovia.com/ggpubr/reference/as_ggplot.html
    1688. #http://www.sthda.com/english/articles/24-ggpubr-publication-ready-plots/82-ggplot2-easy-way-to-change-graphical-parameters/
    1689. #https://plot.ly/r/box-plots/#horizontal-boxplot
    1690. #library("gridExtra")
    1691. #par(mfrow=c(4,1))
    1692. p <- ggboxplot(div_df_melt, x = "Group", y = "value",
    1693. facet.by = "variable",
    1694. scales = "free",
    1695. width = 0.5,
    1696. fill = "gray", legend= "right")
    1697. #ggpar(p, xlab = FALSE, ylab = FALSE)
    1698. lev <- levels(factor(div_df_melt$Group)) # get the variables
    1699. #FITTING4: delete H47(1) in lev
    1700. #lev <- lev[-c(3)]
    1701. # make a pairwise list that we want to compare.
    1702. #my_stat_compare_means
    1703. #https://stackoverflow.com/questions/47839988/indicating-significance-with-ggplot2-in-a-boxplot-with-multiple-groups
    1704. L.pairs <- combn(seq_along(lev), 2, simplify = FALSE, FUN = function(i) lev[i]) #%>% filter(p.signif != "ns")
    1705. my_stat_compare_means <- function (mapping = NULL, data = NULL, method = NULL, paired = FALSE,
    1706. method.args = list(), ref.group = NULL, comparisons = NULL,
    1707. hide.ns = FALSE, label.sep = ", ", label = NULL, label.x.npc = "left",
    1708. label.y.npc = "top", label.x = NULL, label.y = NULL, tip.length = 0.03,
    1709. symnum.args = list(), geom = "text", position = "identity",
    1710. na.rm = FALSE, show.legend = NA, inherit.aes = TRUE, ...)
    1711. {
    1712. if (!is.null(comparisons)) {
    1713. method.info <- ggpubr:::.method_info(method)
    1714. method <- method.info$method
    1715. method.args <- ggpubr:::.add_item(method.args, paired = paired)
    1716. if (method == "wilcox.test")
    1717. method.args$exact <- FALSE
    1718. pms <- list(...)
    1719. size <- ifelse(is.null(pms$size), 0.3, pms$size)
    1720. color <- ifelse(is.null(pms$color), "black", pms$color)
    1721. map_signif_level <- FALSE
    1722. if (is.null(label))
    1723. label <- "p.format"
    1724. if (ggpubr:::.is_p.signif_in_mapping(mapping) | (label %in% "p.signif")) {
    1725. if (ggpubr:::.is_empty(symnum.args)) {
    1726. map_signif_level <- c(`****` = 1e-04, `***` = 0.001,
    1727. `**` = 0.01, `*` = 0.05, ns = 1)
    1728. } else {
    1729. map_signif_level <- symnum.args
    1730. }
    1731. if (hide.ns)
    1732. names(map_signif_level)[5] <- " "
    1733. }
    1734. step_increase <- ifelse(is.null(label.y), 0.12, 0)
    1735. ggsignif::geom_signif(comparisons = comparisons, y_position = label.y,
    1736. test = method, test.args = method.args, step_increase = step_increase,
    1737. size = size, color = color, map_signif_level = map_signif_level,
    1738. tip_length = tip.length, data = data)
    1739. } else {
    1740. mapping <- ggpubr:::.update_mapping(mapping, label)
    1741. layer(stat = StatCompareMeans, data = data, mapping = mapping,
    1742. geom = geom, position = position, show.legend = show.legend,
    1743. inherit.aes = inherit.aes, params = list(label.x.npc = label.x.npc,
    1744. label.y.npc = label.y.npc, label.x = label.x,
    1745. label.y = label.y, label.sep = label.sep, method = method,
    1746. method.args = method.args, paired = paired, ref.group = ref.group,
    1747. symnum.args = symnum.args, hide.ns = hide.ns,
    1748. na.rm = na.rm, ...))
    1749. }
    1750. }
    1751. p2 <- p +
    1752. stat_compare_means(
    1753. method="t.test",
    1754. comparisons = list(c("Group1", "Group2"), c("Group1", "Group3"), c("Group1", "Group4"), c("Group1", "Group6"), c("Group1", "Group8"), c("Group2", "Group5"),c("Group4", "Group5"),c("Group4", "Group6"),c("Group4", "Group7"),c("Group6", "Group7")),
    1755. label = "p.signif",
    1756. symnum.args <- list(cutpoints = c(0, 0.0001, 0.001, 0.01, 0.05, 1), symbols = c("****", "***", "**", "*", "ns"))
    1757. )
    1758. #comparisons = L.pairs,
    1759. #symnum.args <- list(cutpoints = c(0, 0.0001, 0.001, 0.01, 0.05), symbols = c("****", "***", "**", "*")),
    1760. #stat_pvalue_manual
    1761. #print(p2)
    1762. #https://stackoverflow.com/questions/20500706/saving-multiple-ggplots-from-ls-into-one-and-separate-files-in-r
    1763. #FITTING3: mkdir figures
    1764. ggsave("./figures/alpha_diversity_Group.png", device="png", height = 10, width = 12)
    1765. ggsave("./figures/alpha_diversity_Group.svg", device="svg", height = 10, width = 12)
    1766. p3 <- p +
    1767. stat_compare_means(
    1768. method="t.test",
    1769. comparisons = list(c("Group2", "Group4"), c("Group2", "Group6"), c("Group4", "Group8"), c("Group6", "Group8")),
    1770. label = "p.signif",
    1771. symnum.args <- list(cutpoints = c(0, 0.0001, 0.001, 0.01, 0.05, 1), symbols = c("****", "***", "**", "*", "ns")),
    1772. )
    1773. #symnum.args <- list(cutpoints = c(0, 0.0001, 0.001, 0.01, 0.05), symbols = c("****", "***", "**", "*")),
    1774. #stat_pvalue_manual
    1775. #print(p2)
    1776. #https://stackoverflow.com/questions/20500706/saving-multiple-ggplots-from-ls-into-one-and-separate-files-in-r
    1777. #FITTING3: mkdir figures
    1778. ggsave("./figures/alpha_diversity_Group2.png", device="png", height = 10, width = 12)
    1779. ggsave("./figures/alpha_diversity_Group2.svg", device="svg", height = 10, width = 12)
    1780. ```
    1781. # Selected alpha diversity
    1782. ```{r, echo=TRUE, warning=FALSE, fig.cap="Alpha diversity", out.width = '100%', fig.align= "center"}
    1783. knitr::include_graphics("./figures/alpha_diversity_Group2.png")
    1784. selected_alpha_diversities<-read.csv("selected_alpha_diversities.txt",sep="\t")
    1785. knitr::kable(selected_alpha_diversities) %>% kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"))
    1786. ```
    1787. # Beta diversity
    1788. ```{r, echo=TRUE, warning=FALSE, fig.cap="Beta diversity", out.width = '100%', fig.align= "center"}
    1789. #file:///home/jhuang/DATA/Data_Marius_16S/core_diversity_e42369/bdiv_even42369_Group/unweighted_unifrac_boxplots/Group_Stats.txt
    1790. beta_diversity_group_stats<-read.csv("unweighted_unifrac_boxplots_Group_Stats.txt",sep="\t")
    1791. knitr::kable(beta_diversity_group_stats) %>% kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"))
    1792. #NOTE: Run this Phyloseq0.Rmd, then run the code of MicrobiotaProcess.R to manually generate PCoA.png, then run this Phyloseq.Rmd!
    1793. #NOTE: AT_FIRST_DEACTIVATE_THIS_LINE: knitr::include_graphics("./figures/PCoA.png")
    1794. ```
    1795. # Differential abundance analysis
    1796. Differential abundance analysis aims to find the differences in the abundance of each taxa between two groups of samples, assigning a significance value to each comparison.
    1797. ## Group2 vs Group4
    1798. ```{r, echo=TRUE, warning=FALSE}
    1799. library("DESeq2")
    1800. #ALTERNATIVE using ps.ng.tax_most_copied: ps.ng.tax (40594) vs. ps.ng.tax_most_copied (166)
    1801. ps.ng.tax_sel <- ps.ng.tax
    1802. #FITTING5: correct the id of the group members, see FITTING6
    1803. otu_table(ps.ng.tax_sel) <- otu_table(ps.ng.tax)[,c("8","9","10","12","13","14", "21","22","23","24","25","26","27","28")]
    1804. diagdds = phyloseq_to_deseq2(ps.ng.tax_sel, ~Group)
    1805. diagdds$Group <- relevel(diagdds$Group, "Group4")
    1806. diagdds = DESeq(diagdds, test="Wald", fitType="parametric")
    1807. resultsNames(diagdds)
    1808. res = results(diagdds, cooksCutoff = FALSE)
    1809. alpha = 0.05
    1810. sigtab = res[which(res$padj < alpha), ]
    1811. sigtab = cbind(as(sigtab, "data.frame"), as(tax_table(ps.ng.tax_sel)[rownames(sigtab), ], "matrix"))
    1812. sigtab <- sigtab[rownames(sigtab) %in% rownames(tax_table(ps.ng.tax_most_copied)), ]
    1813. kable(sigtab) %>%
    1814. kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"))
    1815. #rownames(sigtab) %in% rownames(tax_table(ps.ng.tax_most_copied))
    1816. #subv %in% v
    1817. ### returns a vector TRUE FALSE
    1818. #is.element(subv, v)
    1819. ### returns a vector TRUE FALSE
    1820. library("ggplot2")
    1821. theme_set(theme_bw())
    1822. scale_fill_discrete <- function(palname = "Set1", ...) {
    1823. scale_fill_brewer(palette = palname, ...)
    1824. }
    1825. x = tapply(sigtab$log2FoldChange, sigtab$Order, function(x) max(x))
    1826. x = sort(x)
    1827. sigtab$Order = factor(as.character(sigtab$Order), levels=names(x))
    1828. x = tapply(sigtab$log2FoldChange, sigtab$Family, function(x) max(x))
    1829. x = sort(x)
    1830. sigtab$Family = factor(as.character(sigtab$Family), levels=names(x))
    1831. ggplot(sigtab, aes(x=log2FoldChange, y=Family, color=Order)) + geom_point(aes(size=padj)) + scale_size_continuous(name="padj",range=c(8,4))+
    1832. theme(axis.text.x = element_text(angle = -25, hjust = 0, vjust=0.5))
    1833. #Error in checkForExperimentalReplicates(object, modelMatrix) :
    1834. # The design matrix has the same number of samples and coefficients to fit,
    1835. # so estimation of dispersion is not possible. Treating samples
    1836. # as replicates was deprecated in v1.20 and no longer supported since v1.22.
    1837. ```
    1838. ## Group2 vs Group6
    1839. ```{r, echo=TRUE, warning=FALSE}
    1840. ps.ng.tax_sel <- ps.ng.tax
    1841. otu_table(ps.ng.tax_sel) <- otu_table(ps.ng.tax)[,c("8","9","10","12","13","14", "33","34","35","36","37","38","39","51")]
    1842. diagdds = phyloseq_to_deseq2(ps.ng.tax_sel, ~Group)
    1843. diagdds$Group <- relevel(diagdds$Group, "Group6")
    1844. diagdds = DESeq(diagdds, test="Wald", fitType="parametric")
    1845. resultsNames(diagdds)
    1846. res = results(diagdds, cooksCutoff = FALSE)
    1847. alpha = 0.05
    1848. sigtab = res[which(res$padj < alpha), ]
    1849. sigtab = cbind(as(sigtab, "data.frame"), as(tax_table(ps.ng.tax_sel)[rownames(sigtab), ], "matrix"))
    1850. sigtab <- sigtab[rownames(sigtab) %in% rownames(tax_table(ps.ng.tax_most_copied)), ]
    1851. kable(sigtab) %>%
    1852. kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"))
    1853. library("ggplot2")
    1854. theme_set(theme_bw())
    1855. scale_fill_discrete <- function(palname = "Set1", ...) {
    1856. scale_fill_brewer(palette = palname, ...)
    1857. }
    1858. x = tapply(sigtab$log2FoldChange, sigtab$Order, function(x) max(x))
    1859. x = sort(x)
    1860. sigtab$Order = factor(as.character(sigtab$Order), levels=names(x))
    1861. x = tapply(sigtab$log2FoldChange, sigtab$Family, function(x) max(x))
    1862. x = sort(x)
    1863. sigtab$Family = factor(as.character(sigtab$Family), levels=names(x))
    1864. ggplot(sigtab, aes(x=log2FoldChange, y=Family, color=Order)) + geom_point(aes(size=padj)) + scale_size_continuous(name="padj",range=c(8,4))+
    1865. theme(axis.text.x = element_text(angle = -25, hjust = 0, vjust=0.5))
    1866. ```
    1867. ## Group4 vs Group8
    1868. ```{r, echo=TRUE, warning=FALSE}
    1869. ps.ng.tax_sel <- ps.ng.tax
    1870. otu_table(ps.ng.tax_sel) <- otu_table(ps.ng.tax)[,c("21","22","23","24","25","26","27","28", "47","48","49","50","52","53","55")]
    1871. diagdds = phyloseq_to_deseq2(ps.ng.tax_sel, ~Group)
    1872. diagdds$Group <- relevel(diagdds$Group, "Group8")
    1873. diagdds = DESeq(diagdds, test="Wald", fitType="parametric")
    1874. resultsNames(diagdds)
    1875. res = results(diagdds, cooksCutoff = FALSE)
    1876. alpha = 0.05
    1877. sigtab = res[which(res$padj < alpha), ]
    1878. sigtab = cbind(as(sigtab, "data.frame"), as(tax_table(ps.ng.tax_sel)[rownames(sigtab), ], "matrix"))
    1879. sigtab <- sigtab[rownames(sigtab) %in% rownames(tax_table(ps.ng.tax_most_copied)), ]
    1880. kable(sigtab) %>%
    1881. kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"))
    1882. library("ggplot2")
    1883. theme_set(theme_bw())
    1884. scale_fill_discrete <- function(palname = "Set1", ...) {
    1885. scale_fill_brewer(palette = palname, ...)
    1886. }
    1887. x = tapply(sigtab$log2FoldChange, sigtab$Order, function(x) max(x))
    1888. x = sort(x)
    1889. sigtab$Order = factor(as.character(sigtab$Order), levels=names(x))
    1890. x = tapply(sigtab$log2FoldChange, sigtab$Family, function(x) max(x))
    1891. x = sort(x)
    1892. sigtab$Family = factor(as.character(sigtab$Family), levels=names(x))
    1893. ggplot(sigtab, aes(x=log2FoldChange, y=Family, color=Order)) + geom_point(aes(size=padj)) + scale_size_continuous(name="padj",range=c(8,4))+
    1894. theme(axis.text.x = element_text(angle = -25, hjust = 0, vjust=0.5))
    1895. ```
    1896. ## Group6 vs Group8
    1897. ```{r, echo=TRUE, warning=FALSE}
    1898. ps.ng.tax_sel <- ps.ng.tax
    1899. otu_table(ps.ng.tax_sel) <- otu_table(ps.ng.tax)[,c("33","34","35","36","37","38","39","51", "47","48","49","50","52","53","55")]
    1900. diagdds = phyloseq_to_deseq2(ps.ng.tax_sel, ~Group)
    1901. diagdds$Group <- relevel(diagdds$Group, "Group8")
    1902. diagdds = DESeq(diagdds, test="Wald", fitType="parametric")
    1903. resultsNames(diagdds)
    1904. res = results(diagdds, cooksCutoff = FALSE)
    1905. alpha = 0.05
    1906. sigtab = res[which(res$padj < alpha), ]
    1907. sigtab = cbind(as(sigtab, "data.frame"), as(tax_table(ps.ng.tax_sel)[rownames(sigtab), ], "matrix"))
    1908. sigtab <- sigtab[rownames(sigtab) %in% rownames(tax_table(ps.ng.tax_most_copied)), ]
    1909. kable(sigtab) %>%
    1910. kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"))
    1911. library("ggplot2")
    1912. theme_set(theme_bw())
    1913. scale_fill_discrete <- function(palname = "Set1", ...) {
    1914. scale_fill_brewer(palette = palname, ...)
    1915. }
    1916. x = tapply(sigtab$log2FoldChange, sigtab$Order, function(x) max(x))
    1917. x = sort(x)
    1918. sigtab$Order = factor(as.character(sigtab$Order), levels=names(x))
    1919. x = tapply(sigtab$log2FoldChange, sigtab$Family, function(x) max(x))
    1920. x = sort(x)
    1921. sigtab$Family = factor(as.character(sigtab$Family), levels=names(x))
    1922. ggplot(sigtab, aes(x=log2FoldChange, y=Family, color=Order)) + geom_point(aes(size=padj)) + scale_size_continuous(name="padj",range=c(8,4))+
    1923. theme(axis.text.x = element_text(angle = -25, hjust = 0, vjust=0.5))
    1924. ```

like unlike

点赞本文的读者

还没有人对此文章表态


本文有评论

没有评论

看文章,发评论,不要沉默


© 2023 XGenes.com Impressum