<jupyter_start><jupyter_text>This notebook is to perform Non-Redundant Analysis on alpha diversities reference: <NAME>, <NAME>, <NAME>, et al. Population-level analysis of gut microbiome variation. Science 2016;352:560-4<jupyter_code>library(vegan)
mf = read.csv('../data/mapping_sleep_alpha.txt', stringsAsFactors=FALSE, sep='\t')
colnames(mf)[1] = 'SampleName'
dim(mf)
head(mf)
#sapply(mf, class)
names(mf)
# correlated variables are not selected, such as the same variable yet one categorical and one continuous
# for example, choose categorical not continuous: 'AMAMPT_C1', 'AMFVT_C1', 'AMPHIT_15SD', 'PQBADSLP', 'PQPEFFIC'

# convert categorial to factors
vars_cat = c('GIERACE', 'SITE', 'GIMSTAT', 'MHDIAB', 'MHRHEU1', 'MHOA', 'MHCHF', 'MHMI', 'MHDEPR', 
             'TUDRAMT_REVISED', 'PQPSLMED', 'PQPSQUAL', 'QLCOMP', 'SLEEPHRS', 'PQBADSLP', 
             'TURSMOKE', 'M1ADEPR', 'M1BENZO', 'PQPEFFIC', 'AMAMPT_C1', 'AMFVT_C1', 'AMPHIT_15SD')
mf[vars_cat] = lapply(mf[vars_cat], factor)

# # convert continuous to numeric
vars_cts = c('Age', 'BMI', 'PASCORE', 'alpha_pd')
mf[vars_cts] = lapply(mf[vars_cts], as.numeric)
length(vars_cat)
length(vars_cts)
dat = mf[c(vars_cat, vars_cts)]
dim(dat)
head(dat)
summary(dat)
dat = dat[complete.cases(dat), ]
print(dim(dat))
alpha = dat$alpha_pd
dat = dat[, -which(names(dat) %in% 'alpha_pd')]
dim(dat)
mod0 <- rda(alpha ~ 1., dat)  # Model with intercept only
mod1 <- rda(alpha ~ ., dat)  # Model with all explanatory variables
step.res <- ordiR2step(mod0, mod1, perm.max = 1000)
table = step.res$anova
table
table$ES.RDA = c(table$R2.adj[1], table$R2.adj[2]-table$R2.adj[1], 
                       table$R2.adj[3]-table$R2.adj[2], table$R2.adj[4]-table$R2.adj[3],
                       table$R2.adj[5]-table$R2.adj[4], table$R2.adj[6]-table$R2.adj[5])
table = table[-6, ]
step.res$call
table
rownames(table) = c('Diabetes', 'Race', 'Site', 'Antidepressant Use', 'BMI')
library(ggplot2)
pdf('../figures/RDA_PD_alpha.pdf')
covariates = rownames(table)
ggplot(table, aes(x=reorder(covariates, ES.RDA), y=ES.RDA, fill=covariates)) +
  labs(x = 'Non-redundant Covariants', y = 'Effect Size (PD Alpha Diversity)') +
  geom_bar(stat='identity') +
 theme(axis.text=element_text(size=10), 
       axis.title=element_text(size=14,face="bold"), 
       legend.position="none") + 
  coord_flip()
dev.off()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.17627261253681112, 'fraction_numerical': 0.02566259991586033, 'mean_word_length': 5.2088772845953, 'pattern_counts': {'":': 0, '<': 11, '<?xml version=': 0, '>': 8, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '1', 'dir': 'jupyter-structured-clean-dedup', 'id': '543886', 'license': 'BSD-3-Clause', 'n_tokens_mistral': 1005, 'n_tokens_neox': 939, 'n_words': 189, 'path': '2.2 Alpha Diversity Analysis -- RDA (R).ipynb', 'repo_name': 'serenejiang/mros_sleep_analysis', 'size': '2391'}