0. Preprocessing.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "library(phyloseq)\n",
    "library(ggplot2)\n",
    "library(vegan)\n",
    "library(readr)\n",
    "library(data.table)\n",
    "library(dplyr)\n",
    "library(stats)\n",
    "library(ggpubr)\n",
    "library(dunn.test)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 79,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "phyloseq-class experiment-level object\n",
       "otu_table()   OTU Table:         [ 729 taxa and 97 samples ]\n",
       "sample_data() Sample Data:       [ 97 samples by 12 sample variables ]\n",
       "tax_table()   Taxonomy Table:    [ 729 taxa by 7 taxonomic ranks ]\n",
       "phy_tree()    Phylogenetic Tree: [ 729 tips and 723 internal nodes ]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "setwd(\"/gpfs/data/gilbert-lab/Lutz/Alzheimers/Update_Oct2019/phyloseq/\")\n",
    "\n",
    "otu <- read.csv(file=\"/gpfs/data/gilbert-lab/Lutz/Alzheimers/Update_Oct2019/phyloseq/final/OTU_X.csv\",header=TRUE)\n",
    "#head(otu)\n",
    "tax <- read.table(file=\"taxonomy.tsv\",sep='\\t',header=TRUE)\n",
    "#head(tax)\n",
    "merged_file <- merge(otu,tax,by.x=c(\"OTUID\"),by.y=c(\"OTUID\"))\n",
    "#head(merged_file)\n",
    "write.table(merged_file, file=\"combined_otu_tax.txt\",sep='\\t',col.names=TRUE,row.names=FALSE)\n",
    "\n",
    "#Export and separate into OTU and taxonomy tables (see qiime2_to_phyloseq.pdf)\n",
    "\n",
    "otu_table=read.csv(\"OTU_table.csv\",sep=\",\",row.names=1)\n",
    "otu_table = as.matrix(otu_table)\n",
    "taxonomy= read.csv(\"tax_table.csv\",sep=\",\",row.names=1)\n",
    "taxonomy=as.matrix(taxonomy)\n",
    "\n",
    "OTU=otu_table(otu_table,taxa_are_rows=TRUE)\n",
    "TAX=tax_table(taxonomy)\n",
    "\n",
    "OTU=otu_table(otu_table,taxa_are_rows=TRUE)\n",
    "TAX=tax_table(taxonomy)\n",
    "\n",
    "META=import_qiime_sample_data(\"/gpfs/data/gilbert-lab/Lutz/Alzheimers/Update_Oct2019/phyloseq/final/META_X.txt\")\n",
    "tree=read_tree(\"rep_phylo_tree.nwk\")\n",
    "\n",
    "ps = merge_phyloseq(OTU,TAX,META,tree)\n",
    "ps\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 80,
   "metadata": {},
   "outputs": [],
   "source": [
    "alz_nt <- prune_samples(sample_sums(ps) > 1000, ps)\n",
    "alz_t <- merge_phyloseq(OTU,TAX,META,tree) %>%\n",
    "    prune_samples((sample_sums(.) > 1000), .) %>%\n",
    "    filter_taxa(function(x) sum(x > 3) > (0.001*length(x)), TRUE) %>%\n",
    "    transform_sample_counts(function(x) x / sum(x))\n",
    "\n",
    "saveRDS(alz_nt, \"/gpfs/data/gilbert-lab/Lutz/Alzheimers/Update_Oct2019/PS_for_analysis/alz_nt.rds\")\n",
    "saveRDS(alz_t, \"/gpfs/data/gilbert-lab/Lutz/Alzheimers/Update_Oct2019/PS_for_analysis/alz_t.rds\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Warning message in strsplit(msg, \"\\n\"):\n",
      "“input string 1 is invalid in this locale”"
     ]
    }
   ],
   "source": [
    "#setwd(\"/gpfs/data/gilbert-lab/Lutz/Alzheimers/ps_files/\")\n",
    "\n",
    "#biom <- import_biom(\"abx_deblur.biom\")\n",
    "#map <- import_qiime_sample_data(\"ABX_mapfile.txt\")\n",
    "#tree <- read_tree(\"rep_phylo.tre\")\n",
    "#alz_nt <- merge_phyloseq(biom,tree,map)\n",
    "#alz_nt <- prune_samples(sample_sums(alz_nt) > 1000, alz_nt)\n",
    "#alz <- merge_phyloseq(biom,tree,map) %>%\n",
    "#    prune_samples((sample_sums(.) > 1000), .) %>%\n",
    "#    filter_taxa(function(x) sum(x > 3) > (0.001*length(x)), TRUE) %>%\n",
    "#    transform_sample_counts(function(x) x / sum(x))\n",
    "\n",
    "#saveRDS(alz_nt, \"/gpfs/data/gilbert-lab/Lutz/Alzheimers/RDSfiles/alz_nt.rds\")\n",
    "#saveRDS(alz, \"/gpfs/data/gilbert-lab/Lutz/Alzheimers/RDSfiles/alz.rds\")                            "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 73,
   "metadata": {},
   "outputs": [],
   "source": [
    "#alz_nt = readRDS(\"/gpfs/data/gilbert-lab/Lutz/Alzheimers/RDSfiles/alz_nt.rds\")\n",
    "#alz = readRDS(\"/gpfs/data/gilbert-lab/Lutz/Alzheimers/RDSfiles/alz.rds\")\n",
    "\n",
    "#Read in NEW DATA NOVEMBER 2019\n",
    "alz_nt = readRDS(\"/gpfs/data/gilbert-lab/Lutz/Alzheimers/Update_Oct2019/PS_for_analysis/alz_nt.rds\")\n",
    "alz_t = readRDS(\"/gpfs/data/gilbert-lab/Lutz/Alzheimers/Update_Oct2019/PS_for_analysis/alz_t.rds\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 81,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "`stat_bin()` using `bins = 30`. Pick better value with `binwidth`.\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA0gAAANICAIAAAByhViMAAAABmJLR0QA/wD/AP+gvaeTAAAg\nAElEQVR4nO3daZxbZaHA4fckmWTW7mXfKWVRZBUKAlKpCCi77LIvIkJFtgqyKJsFWYR7AYWC\ngFwQRfQiCl5EQFEWBUVBQEBR9qXQZaazZnI/pAyly3TKpEnmzfN86C9zcnLOe3KSk39PMpmk\nUCgEAACGvlSlBwAAQGkIOwCASAg7AIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCASFR12BV6\n22+/8hs7bbPx2JEtmVRm2OiVJnx6zwtv/HVP7X2ncpIkSZKUf43zSmWyo1dcc4d9vnTLb/+5\nVFf94IMPPvjggwsdz1JdLwAMdUnV/uWJ3u43v7T1x65+5I0kyay+3gYrjsxNf/m5v7/4Vghh\ntU+f/Le7LmhO19DLfLFpyrmzimtce+21+6Z0d8x++aXXunoLIYRJX7nuV985dCn9t2ChG1v+\newAAhpzqDbuHT9loi2//ZdT6B/78jiu2XK2lOPGlR3+y/04HPji9faupf/ndlA0qO8Jyam1t\nDSE0NzeXbY0LDanu1lfuuO6SY06+7M2u/KfOfejer08o26qFHQAsVvWG3fajGu55t+Mnb87Z\nY2zDvNOn//X0MRuc1zj2821v/rhSY6sF/YTUmw9fstInTsqnhv3+nTcntGTLs2phBwCLVb2f\nsXtwZmcIYbtR9fNNH/WRKZtsssl6q78+3/QZz/76uAM+u/ryY7LpzMjl19hhv2PvePzNBRf7\nyK0X7Lj52sMb65pGr7LX5Eve7umd78Nbi/os10KnD2Sl792w9/ZvHb3OiiPS6ezy4zY+buqt\n3QskyisP/+jQ3bZdYeyITDq73BofO/CEC59r7V7U2ge+2L/89JLPbbnuyKa6ptGr7nviFW29\nhUF+Xm2ZCSdcv/1KvT0zJ1/17BLdG++tt/CbaWd9coM16uvSjcOW2Xz7fa+957l551nUJocQ\nCvnWK07ab81lhqXTdcuuscEXz76xU+kBQJ9CtdpyWC6EsOt37hvIzP++8+zhmVQIYfnxG227\n7dbjl28KISSp3NHfe3ze2e458zPFrV71ox+fsPH4dJKssO1Z890Pi7pbFpw+wJUWb/jQeZMy\n9SvssMfeO2/7seKUzb/2u3lne/LaY+pSSQhhxfU2/eQ2W6zYXBdCaFhmq8dmdy107QNc7INT\n9wghJEl67U0+sdXH16tLkrX2nTaQXd//PG//7UshhOGrnbdE90ZxmT/76pYhhJaV1vnEVlus\nMbY+hJAkyX4XP1ScZ8qUKcXZpkyZMmXKlPlu+53d18g0rPy5vfff5VMbFLNv4+N/0/+GAEDt\nqN6w++PZnyy+lq+0yWdOveC7v3vixfwi5uyc9fAq9Zm6xvWu/+2/5k7q7f7tjacPy6SSVO7a\nF2cVp83+zzV1SZLKDL/orn8Up7z6h2uXy6Y/XNgNcKV9N1xmzCcfeau9OOVPV342hFDX9NG+\neea89b/DM6lUZvhFv3ymOKWn47Uzt18phLDWgb9e6KgGstjW1/6nPpWks8tf99CrxSlvPn7z\nqvWZwYdd56xHQgjZ5g2X6N4oLjNJ1R0/7f65ezPf/r8XHZwkSZLK3fRKaz+rLk6sH7ntn97p\nKE554oZ9Qgh1jev2vyEAUDuqN+wKhZ4fnf/FNUfmwnuall9794MnX3fHI/MV3qMnrR9COOSe\nl+e7/eMXbRlCWGnS7cUf79pt9RDCBqf+Yd55Hjt3kw8XdgNcad8N9/jFv/um5DtfnW9p9x24\n1oJj65hxfwihYfSuCx3VQBZ7z95rhhA2u+DP897wr5dMGHzY9ebnhBCSVH3xxwHeG8VlrnPU\nL+eb7dZ91gwhrL7bXf2sujjx83f95/3t7Z4+kA0BgNpR9S+KvZ1P/uGui886fvsJ62ZTcz9x\ntcrEY17tfL/ujl6+OYTwWtf8Z/Ta3749hNAwaqfij7uPaQgh3PRG27zzzHnzlg8XdgNcad8N\nH5/d1c/S9hnbGEL4wQfHVigU/vPiiy+++PJCRzWQxX5+bGMI4ba35nxgk9+6bfBhV5gbdg3F\nnwZ4bxSXOe31+Tdz9suXhRAaRu/cz6rnbm9rf9sLADWuen8rdkGd77740xuumHLqd/7T0bP2\nYXc/c+3cD8ytUl/3UmfPom6Vrhvb0/VmCGGZbOat7vysnt6Web4Ar5CflcoMD/P8uuWifvty\nvukDXGnfDfOFwry/qDLf0pbLZd7oys/o6R2+iC/nW3BUA1lscZNb871NqXk3eXYqM2yh29j/\nGufV3frHbMtm2eaNOmc/HgZ8bxSXOd8uCCEUet5J1Y1OZYblu2cuatUD2V4AqHGZSg9g4a69\n9tpCoXDYEUfM+yqeG7navsd/e7uthy2z6Zn/vPXE8F7YvdGdDx/8Kt15JcncT9HNyveGEPLz\nNcB71y5Ofr6fB7jSPv3/+vHMnt7FzrNQ/d+kuMnz/6GOgW5yv0t+8foQQuPYvYs/LtG9seDv\n7YakLoRQ6O1e7Hqr97e4AaAKVGnYfePLR7/c2bPO3gdtNWz+r0kbs9GJIZzZ0/7+X7Vav6nu\nsdldv/vb38fW9fe6v0FT3aOzu37+TvuByzT2Tex45+6BjKdr9h/nmzLAlQ7Q5i3ZB2Z2/npG\nx+6jP/Clfb+55ab/dPbsddDBH26xH2/JPjiz89fvduw55v3Fds64b1BjDSGE8H+n3BlCGP/l\nzxV/XKJ747a35hy1fNO8U4pviOeGbTH4gQFALavSMyDHrj4shHD013+54FWz/nlFCKFxzB59\nUyZ/ZFQI4fTfzf/Ndv/+3+MnTZq0z1fm/tXRo8ePCCFcdOlf5p3n2WkXLHQAL3d94BTdv358\nxnwzDHClA3TcpmNDCOdf9dS8E3vmPLnzgQcfecy5DakP+Z1zx24wOoRw4XUf+La5f9268E0e\nuDcfvuTgX72Uyoz4r6PGF6cs0b1xyZn3zzfb/025KISwzITJgxwYANS6yny0b3He/vPUTJKE\nED795Qv//OK7fdPfeOqeXVZpCSFMvPzJvonvPDU1SZK6xrWvve+5vonP3fu91eozIYRTHn5j\n7m0fOSGE8IGvO3n4ugW/7mTdxroQwjZnv/+bm68+csOq9bn5ZhvgSgv9/ipA348zX7ginSSp\nulEXv/d1J/mu6efvtmoIYc2971joQgay2HeeviCEkM6t0Pd1J9P/+uM1GuoGsusXOk9322u3\nXX7Sstl0CGHS+Y8s6b1RXGaSqpvc93UnvZ13XnpYKkmSVN21//7At6LMyfcu6fYCQI2r3hfF\nx6/76sjM3BOKw5dZaa3xa6+6/Mjijyt88oS2D77q33XmZ4tXLTd+o20/NXGj8csVf/zUV2+d\nd7Yr9l2nOH2Vj266xSbrpJNkmc1Pni8O/nDGNnPX8rFt9tp33x22Xj+dJDufe/+CDTHAlQ6w\nSO47b/e5611n42233Wq1EbkQQuMy2zzR2t8XFC92sT87cZsQQpKkx2+8xTYT1s+mko9Pvj2E\nkKRyi77v31/U2vMYt/pKuffOHU76ynXz/QbsQO6N4pSvfnL5EELLSut+YutPjFumIYSQJMk+\n3/5932wr5NIhhHW33/2A/fZe0u0FgFpW1S+Ks//98NSTj9hyg/EjmrIhhFzzqPW32OHMK346\nZ2FfVfzUL7637/abLzOsPkmSltErbrXTAVf/4m/zz9Tb9ZOLT9h8nRWzqSQ3fPkdD/vGfzp6\nFoiD/F1Xnrb1x9ZoyKRS6fo1N/rUt27+Y6FQSKfT6XT6Q6x04EXy9198b6/tNh3dnEuS9JiV\n193n2HOfntW1qPkHvtjf3XDexA1Xz6aS5mXWPPzsWzq7pocQMg1rzH/nLGxR80qS9PBlV5u0\n51E33ff8Qm+y2HujuJzentnXn3PMRmsuV5dKsk2jN9lu72t+9Y95Z3vimhPWXm5YCKFl9Iof\nYnsBoGYNpa87WUoi/sqMzrdff6ujp3G5FUZl3v8wZcf0/20Ys1vTsge2vn5jmccT8V0NANWg\nSn95gpL424WfWXnllT/73Wfmnfj8LReFEEZvuE+FBgUALC3CLmZrH3tECOFPJ+908W0PzujI\n97S9fc+NZ2934kMhhC9eumWlRwcAlJi3YiN/f/DHk7fa+79+P9/ErY774e8ur8AZu7jvagCo\nuCr9gmJKZa/LH3zkk1f89w9+9sRzr+Qzjauvt+keBx1z6I4fq8hgrrnmmoqsFwBqhDN2AACR\n8Bk7AIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCASAg7AIBICDsAgEhU3V+eKBQKs2bN6n+e\nXC7X3d3d29tbniFVibq6unQ63dvb29XVVemxlFWSJLlcrqOjo9IDKbdcLpckSU9PT09PT6XH\nUlbpdDqTyXR2dlZ6IOVWX18fQujq6qrBg1sIobu7u9IDKatUKpXNZkMItXlw6+npyefzlR5I\nWWUymUwmUygUSnJwGz58+CJXNPill9xin94tLS1dXV21dhTI5XLFom1ra6v0WMoqlUoNGzas\ntbW11v5KSlNTUyaT6e3trbWHeiqVymQyra2tlR5IWRX/AxNC6OzsrLU9XtzwWtvqbDZb3PDa\nPLh1d3fX5h7v6elZ2gc3b8UCAERC2AEARELYAQBEQtgBAERC2AEARELYAQBEQtgBAERC2AEA\nRELYAQBEQtgBAERC2AEARELYAQBEQtgBAERC2AEARELYAQBEQtgBAERC2AEARELYAQBEQtgB\nAERC2AEARELYAQBEQtgBAERC2AEARELYAQBEQtgBAERC2AEARELYAQBEQtgBAERC2AEARELY\nAQBEQtgBAERC2AEARCJT6QEAVWHy5MklWc7ll19ekuUA8CE4YwcAEAlhBwAQCWEHABAJYQcA\nEAlhBwAQCWEHABAJYQcAEAlhBwAQCWEHABAJYQcAEAlhBwAQCWEHABAJYQcAEAlhBwAQCWEH\nABAJYQcAEAlhBwAQCWEHABAJYQcAEAlhBwAQCWEHABAJYQcAEAlhBwAQCWEHABAJYQcAEAlh\nBwAQCWEHABAJYQcAEAlhBwAQiUylB7AQI0aM6H+GVCrV1NTU2NhYnvFUiVQqFULIZDKLvX+i\nNHz48EoPodzS6XQIIZfL1dXVVXosS2Dwj88kSVKpVG0+zkMINXtwq7U9niRJ8UJtHtwaGxsb\nGhoqPZCyKj7O0+n04B/qvb29/VxbjWHX0dHR/wxNTU3d3d09PT3lGU+VyOVyqVSqt7d3sfdP\nZJIkyWQynZ2dhUKh0mMpq8bGxiRJenp6urq6Kj2WJTD4x2cmk8lms7X5OA8h1ObBLYTQ2dlZ\n6YGUVSaTKf7nrTYPbjX4OM9ms9lstlAoDP7gVigUstnsoq4dkmHX2NjY1dU1tF7tBi+TydTV\n1dVg2BVP0HZ0dNTasa++vj6VSuXz+aG1xwc/2lwul8lkhtZWD16SJE1NTSGE2jy4hVI8coaW\nbDZbX18fQqjNg1t3d3etpXyx40v1It7S0rKoq3zGDgAgEsIOACASwg4AIBLCDgAgEsIOACAS\nwg4AIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOACASwg4AIBLCDgAg\nEsIOACASwg4AIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOACASwg4A\nIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOACASwg4AIBLCDgAgEsIO\nACASwg4AIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOACASwg4AIBLC\nDgAgEsIOACASwg4AIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOACAS\nwg4AIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOACASwg4AIBLCDgAg\nEsIOACASwg4AIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOACASwg4A\nIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOACASmfKsppBv++nVl/zvA0/M7m3cZNI+U478bCYp\nz5oBAGpFmc7YPfHdU376wvBjz5g69Ywv1T918zfveqk86wUAqB1lOmM37YFXT73+svUaMyGM\nO+683BeOvTPs9KW+a99555329va+HxsbGxe7wFQqlU6nl8pYq1WSJMV/a3PD0+l0oVCo9FjK\naoju8cGPNpVKDbmtHrzi7g41fHCrta1OpeaeWKnNg1vNPs5LcnDr/wFTprB7u7t3bN3cB3Em\nt3LXrO+E8H7YXXLJJXfffXfxcnNz8/3337/YBTY3Ny+FYQ4BmUxm5MiRlR5FBYwYMaLSQ6iM\n+vr6+vr6So9iCZTq8Vmbj/NQwwe3XC5X6SFURm0e3Jqampqamio9igpIp9ODP7j19vb2c22Z\nwm7PlVrOv+6eMw/6VFPPO3d9/9xCfmZ51gvUsv33378ky7n55ptLshyApa1MYbfrOSf/9dSp\nB//iiiSpm/D5ySH1vXmvPeaYYw444IDi5VQqNWPGjP6XNmzYsPb29u7u7qU13KrU0NCQy+V6\nenpaW1srPZaySqVSw4YNmzlzZq29W9HS0pJOpzs6Ojo6Oio9liWw2OfvYmWz2VwuN3v27JKM\npyQGv1GLlSTJ8OHDQwhtbW21dnArfvxmzpw5lR5IWdXV1RVPWdXmwa2zs7Orq6vSAymr4tsv\n+Xy+JAe3fk70linssiM2OOeqm99+/dXQtMzw1N/3u3udea9dYYUVVlhhheLlQqEwffr0xS4w\nn8/39PQslbFWq+Izv1Ao1NqGFz+J0tPTU2vHviG6xwc/2uJHjqpqq8swmL7P2NXgwa34vlKt\nbXXfZ+xq8+BWg4/zsh3Sy/RbsQ/890V3vdUxZrkVx7TUvf3Hnwxbc/vyrBcAoHaUKexWHvv6\njWf+95+e+sffHvrF2Vf+fZ9jP1ae9QIA1I4yvRW7+l5n7/yfsy88/eTUiFV2/vLFn1lm8V9o\nAgDAEilT2CWpxv1Pnlqa308DAGBh/K1YAIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCASAg7\nAIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCASAg7AIBICDsAgEgI\nOwCASAg7AIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCASAg7AIBI\nCDsAgEgIOwCASAg7AIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCA\nSAg7AIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCASAg7AIBICDsA\ngEgIOwCASAg7AIBICDsAgEgIOwCASGQqPQCA+R100EGVHgLAkOSMHQBAJIQdAEAkhB0AQCSE\nHQBAJIQdAEAkhB0AQCSEHQBAJIQdAEAkhB0AQCSEHQBAJIQdAEAkhB0AQCSEHQBAJIQdAEAk\nhB0AQCSEHQBAJIQdAEAkhB0AQCSEHQBAJIQdAEAkhB0AQCSEHQBAJIQdAEAkhB0AQCSEHQBA\nJIQdAEAkhB0AQCSEHQBAJIQdAEAkMpUewELU19f3P0OSJNlsNpWqrSpNp9MhhFQqtdj7JzJJ\nkoQQ6uvrC4VCpcdSVsVHeDqdHlp7fPCjzWSq7rhUhl1QfJyHEGr24Da0HueD1/c4r82DW11d\nXd9jvkaU8EW8/wdM1R1Aw8DCrq6urgqP/ktV8Vhfg2FXlMvlKj2Ecise9TKZzNB6mR/847MK\nD/flfNLV7MGt1o5sfY/z2jy41ezjPEmSwT/Ue3t7+7m2Gu/WGTNm9D/DqFGj2traurq6yjOe\nKtHc3FxfX9/T0zNz5sxKj6WsUqnUqFGjZs6cWWv/qR0xYkQmk+ns7Gxra6v0WJbAYp+/i1WF\nr3OD36jFSpJk9OjRIYTaPLiFEFpbWys9kLLKZrPDhg0LIdTmwa29vb2zs7PSAymrpqamhoaG\nfD5fkuPJmDFjFnXVUDoTAABAP4QdAEAkhB0AQCSEHQBAJIQdAEAkhB0AQCSEHQBAJIQdAEAk\nhB0AQCSEHQBAJIQdAEAkhB0AQCSEHQBAJIQdAEAkhB0AQCSEHQBAJIQdAEAkhB0AQCSEHQBA\nJIQdAEAkhB0AQCSEHQBAJIQdAEAkhB0AQCSEHQBAJIQdAEAkhB0AQCSEHQBAJDKVHgCUz+TJ\nk0u1qMsvv3yQS6iqwZRQCbeLhaq2e7iqHn6AM3YAAJEQdgAAkRB2AACREHYAAJEQdgAAkRB2\nAACREHYAAJEQdgAAkRB2AACREHYAAJEQdgAAkRB2AACREHYAAJEQdgAAkRB2AACREHYAAJEQ\ndgAAkRB2AACREHYAAJEQdgAAkRB2AACREHYAAJEQdgAAkRB2AACREHYAAJEQdgAAkRB2AACR\nEHYAAJEQdgAAkRB2AACREHYAAJEQdgAAkRB2AACREHYAAJEQdgAAkRB2AACREHYAAJEQdgAA\nkRB2AACREHYAAJEQdgAAkRB2AACREHYAAJEQdgAAkRB2AACREHYAAJEQdgAAkRB2AACREHYA\nAJEQdgAAkRB2AACREHYAAJEQdgAAkRB2AACREHYAAJEQdgAAkRB2AACRyJRtTU/8/Jqrb7v3\n5Vm9a2366VNOOXyZOk0JAFBKZaqrd/76vXNu/ecex379ovNPXafj92dd8dfyrBcAoHaU6Yzd\nfVc+sNPZ391ujWEhhDVOP2eDp/LlWS8AQO0oU9j9+u32I5+7/cvfvPuV2eEjn9h5yvH7z3vt\n97///T/+8Y/Fyw0NDVOnTu1/aUmSNDY2NjQ0LK3hVqV0Oh1CyGQyw4cPr/RYKmDYsGGVHsIH\nVNVeqKrBRKmc9/CQO7gN/s4pHtxq7WGcJEnxQrUd3MognU43NjbW19dXeiBlVXycp9PpwT/U\ne3t7+7m2TGH3WnfvpXe8cczJ5y6bmXXn1Zec+r2PXnHMBn3XvvDCC48++mjxcnNzc11d3WIX\nmMmU79OBVSVJkoHcP/Gptq2uqvFU1WCidPDBB5dkOTfffPNi5xlyB7dSPfxSqRr94HVtPn/T\n6XQxdGpNSV7EqyLsMkk44tzjtxiZCyF88awv7nfk9eGYS/uu/fjHP97Y2Fi8nMvlOjo6+l9a\nLpfr7u7uf8PiU1dXl06ne3t7u7q6Kj2WskqSZCCPijKrqvFU1WDoR/97qngCo6ura2gd3Ab/\n8Cu+znV3d5diOENGKpXKZrOhJp+/2Ww2n8/n87X1oaxMJpPJZAqFQmdn5yAXVSgU+jmvX6aw\nWy2XnjAiW7ycG7Z5vvPSea/dddddd9111+LlQqEwffr0/peWzWY7OjpqrW+am5vT6XQ+n29t\nba30WMoqlUrlcrm2trZCoVDpsbyvqvZCVQ2GfvSzp5IkKYbdkDu4Df7h19zcXJLlDC3ZbLYY\ndtV2cCuDESNGdHR0DL5vhpampqZMJlOqF/F+wq5Mp763W6bxF6/NKV5uf/veuuYNy7NeAIDa\nUaYzdpsf9YmjT5066iv7LZ9+98f/df34PRfz6xEAACypMoXdqPW/dNynL5o29fTZoWWzzxx6\n9u6rl2e9AAC1o3y/frX1F07a+gtlWxsAQM2p0V8vBwCIj7ADAIjEQMMuSZIvPz9jwekv3bN9\nXeOaJR0SAAAfxmDP2L315Mx8x79LMhQAAAZjMWF35VqjkiQp/km7K9camSxgkxMerWtavyxD\nBQCgP4M9Y5ekcjt97bslGQoAAIOxmK87Oea5d44JIYSQJMkxz717xbgRZRgTAAAfwkDP2OW7\nXms4ZcfmbHrBd2OLb9QCAFBZA/2C4geP3+7in/59qQ4FAIDBGOgZu7Nu/mdu2Ba//OvL+cJC\nLNUhAgAwEAM9Y/fQrM4NL7lmx/VXXKqjAQDgQxvoGbvxDZmGFRuW6lAAABiMgYbdJQeNf+pb\nt3jPFQCgag007La7/MGjVrl5wqHnPPzMy115gQcAUHUG+hm7VN3IEEIIZ25x/ZkLXuv3JwAA\nKm6wf3kCAIAqMdAzds7JAQBUuYGGXf9/XkL2AQBUnLdiAQAiMdCwm/9vTeS7Xn32kXMP/8SK\n2xzzWld+qQ4RAICB+LBn7FJ1y4/f7OvTHjx2zm2fPPLukg4JAIAPY7Bvxe5z2kf+9eOTSjIU\nAAAGY7Bh1/5ye0/H8yUZCgAAgzGosHvzmV8fdcafsy0TSjUaAAA+tBJ83ck6R3yjJEMBAGAw\nBhp2C5WkchvscOTPr/hkqUYDAMCH5i9PAABEYgk+Y5fv+M/lpx854aOrD2vMpesaV15708O+\n9p0X23uW3uAAABi4gZ6x62578rNrT7jnlba+KS//47HvX/DYj26685F/3P2RxkG9pQsAwOAN\n9Izdb47e455X2rY45MzfPP6P2R3dPZ2tzz3+mzMO2rztlXt3O/LXS3WIAAAMxEDPtJ19x3+W\n3eyiP3z/xL4bjtto4tk3PNz892XP+PmZIeywlMYHAMAADfSM3eOtXavu/ZkFp++032rdbX8p\n6ZAAAPgwBhp26zXWvXjLLxec/qtb/13XtEFJhwQAwIcx0LD75k4rvfnYlAkHn3Xfn59v6+zp\n7Z7zzyd++83DPnHSo2+svOM3luYIAQAYkIF+xu4z194+8XefuO/Gsz9149nzTm9cbuJPrt1+\nKQwMAIAlM9AzdnXNG/zfc09ccsqhG6+1Yi6dSqVzK6610SEnX/zX5/9vg+a6pTpEAAAGYom+\noHhOdviye1/3WEdPPt/Tcftn1qtP97Z15pfe4AAAGLiBhl37m3dtuOomx3596i//Nas4pe2l\nP3136skfX3WTu95sX2rDAwBgoAYadrd+/ojnukdd8au/3XfgWsUpE3/2zHO/uXp01zOH73HL\nUhseAAADNdCwO+dPb66x9y3HbP/ReW8wbuKRt+y/5luPnbc0RgYAwBIZaFZLW/IAABykSURB\nVNi90pkf/pHRC04ftd7wfNcrJR0SAAAfxkDDboeR9c9Nu76nMN/k3pumPV8/YrtSjwoAgCU2\n0LA787j1Z75w6fp7nHDXI0/P7uju7en4119/e/r+m1z4j3fHH3bGUh0iAAADMdAvKN7o63cf\nffc63/3ZpTv97NJ5p49cd5+7zttsKQwMAIAlM9Azdklm5JUPPn/LRSdts8GaTdl0Kp1bcfwm\nXzzzqqf/cvPy2SX4MjwAAJaSgZ6xCyEk6ZZ9T/z2vid+e+mNBgCAD83JNgCASAg7AIBICDsA\ngEgIOwCASAg7AIBICDsAgEgIOwCASAg7AIBILMEXFANANZs8eXJJlnP55ZeXZDlQfs7YAQBE\nQtgBAERC2AEARELYAQBEQtgBAERC2AEARELYAQBEQtgBAERC2AEARELYAQBEQtgBAERC2AEA\nRELYAQBEQtgBAERC2AEARELYAQBEQtgBAERC2AEARELYAQBEQtgBAERC2AEARELYAQBEQtgB\nAERC2AEARELYAQBEQtgBAERC2AEARELYAQBEQtgBAERC2AEARCJT6QEsxOjRo/ufIUmSlpaW\n8gymeiRJEkKoq6tb7P0TpVGjRlV6CB9QVXuhqgZDPwayp4bcwW3wD7/iwS2Xy5ViOKVRzudU\ntR3cyiBJkubm5ubm5koPpKyKj/NMJjP4R1c+n+/n2moMu1mzZvU/Q0tLS3t7e09PT3nGUyUa\nGhqy2WxPT09bW1ulx1JWqVSqpaVl9uzZhUKh0mN532IfpeVUVYOhH/3sqSRJhg0bFkKYM2fO\n0Dq4Df7h19DQEEJob28vxXBKowzPqUwm09TUFEKotoNbGTQ3N3d2dnZ3d1d6IGVVX1+fy+Xy\n+Xxra+vglzZ8+PBFXVWNYTeQnZ3P52vtMVH872yhUKi1DU+lUiGE7u7uqjr2VdVeqKrB0I9+\n9lTxf/NhCB7cBj/a4sGtqra6DIPp2+PVdnArg0KhMOQe54OXzWZDWV7EfcYOACASwg4AIBLC\nDgAgEsIOACASwg4AIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOACAS\nwg4AIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOACASwg4AIBLCDgAg\nEsIOACASwg4AIBLCDgAgEsIOACASmUoPABisyZMnV3oIAFQFZ+wAACIh7AAAIiHsAAAiIewA\nACIh7AAAIiHsAAAiIewAACIh7AAAIiHsAAAiIewAACIh7AAAIiHsAAAiIewAACIh7AAAIiHs\nAAAiIewAACIh7AAAIiHsAAAiIewAACIh7AAAIiHsAAAiIewAACIh7AAAIiHsAAAiIewAACIh\n7AAAIiHsAAAiIewAACIh7AAAIiHsAAAiIewAACIh7AAAIiHsAAAiIewAACIh7AAAIiHsAAAi\nIewAACIh7AAAIiHsAAAiIewAACIh7AAAIiHsAAAiIewAACIh7AAAIiHsAAAiIewAACIh7AAA\nIiHsAAAiIewAACIh7AAAIiHsAAAiIewAACIh7AAAIiHsAAAiIewAACIh7AAAIpEpz2raXnr4\niu/+8M/P/qe9N7PKupsddsKxG46uL8+qAQBqRDnO2BXyrWecdOGzufVP/ObUSy4461MrvzH1\nazeWYb0AADWlHGfsOmfc+2qSPfOkg9drzIQQVl/jzO/fdWgIR5Vh1QAAtaMcYVc/etcf/nDX\nEELHm6+9Pnvms4/e3rTijvPO8OSTT77++uvFy+l0etNNN+1/gUmS1NXVJUmylAZcndLpdAgh\nlUrlcrlKj6Wsijs6m81WeiAfUGt7gZLo52HTd0Abcge3wT8Xige3qnpOlWEwmUymb12FQmFp\nr66qJEnSt/m1o4Qv4v0/YMp6zz757a+d/ey7SZI+8OIb5p3+wx/+8O677y5ebm5uvv/++xe7\nqIaGhqUxwuqXTqdbWloqPYoKqLatrrbxMCR88YtfrPQQSq9Uz4W6urqSLKckyvkEb25uLtu6\nqkfNvoinUqnBP7p6e3v7W8Ugl75ENv32DT/98f+csv9Hf3TODYufGwCAJVGOM3ZzXnnyydey\nm206PoSQzrVs+flTLrz50BCO65vh61//+imnnNL34/Tp0/tf4MiRI9va2rq6upbSgKtTU1NT\nfX19T0/PzJkzKz2WskqlUiNHjnznnXeq6t2KxT5KoUYM/rlQPGXV2tpaiuGURhme4Nlstnjm\nptoObmUwfPjwjo6Ozs7OSg+krBobGxsaGvL5/IwZMwa/tNGjRy/qqrKE3Rs/Of/8F2748Q3D\n00kIId/1SqruAwNqaGjoOytbKBQG8owqFAq19kwoqsENL25vtW14VQ0GKmjwz4W+53gphlMa\nZRhM3yqq7eBWHrW51aEsG16Ot2JHfuSwManZJ19ww2NP/eO5px///vlTx252SBnWCwBQU8px\nxi6dW/mSqV+5YtptF5zxs65Uw3qb7/itr04ow3oBAGpKmX4rdti4bU+dum151gUAUJv8rVgA\ngEgIOwCASAg7AIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCASAg7\nAIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCASAg7AIBICDsAgEgI\nOwCASAg7AIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCASAg7AIBI\nCDsAgEgIOwCASAg7AIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCA\nSAg7AIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCASAg7AIBICDsAgEgIOwCASAg7AIBICDsA\ngEhkKj0ASmzy5MmlWtTll19eqkXFp4T3M3DYYYdVegil52hMRThjBwAQCWEHABAJYQcAEAlh\nBwAQCWEHABAJYQcAEAlhBwAQCWEHABAJYQcAEAlhBwAQCWEHABAJYQcAEAlhBwAQCWEHABAJ\nYQcAEAlhBwAQCWEHABAJYQcAEAlhBwAQCWEHABAJYQcAEAlhBwAQCWEHABAJYQcAEAlhBwAQ\nCWEHABAJYQcAEAlhBwAQCWEHABCJTKUHsBBNTU39z5AkSX19fV1dXXnGUyUymUwIIZ1OL/b+\nKZWyrah/SZKEEJqamgqFQqXHAsxv8AeK4sGtqlTJ0a9PtY1nkFKpVC6Xq8L9vlQVoyWVSg1+\nb/b29vZzbTXeranU4s8jJkkykNliUuybMLD7pySq5B4ubniSJH33AFA9Bn+gqMKndpUc/fpU\n23gGqXg8j2yjFqvvtWxpP2WqMexmz57d/wyjRo1qb2/v6uoqz3iqRHNzczqdzufzi71/SqVs\nK+pfKpUaNWpUa2urM3ZQhQZ/oGhubi7JSEqoSo5+faptPIM0YsSIjo6Ozs7OSg+krJqamhoa\nGkr1Il5fX7+oq2qrlwEAIibsAAAiIewAACIh7AAAIiHsAAAiIewAACIh7AAAIiHsAAAiIewA\nACIh7AAAIiHsAAAiIewAACIh7AAAIiHsAAAiIewAACIh7AAAIiHsAAAiIewAACIh7AAAIiHs\nAAAiIewAACIh7AAAIiHsAAAiIewAACIh7AAAIiHsAAAiIewAACIh7AAAIpGp9AAqafLkyZUe\nwgdcfvnllR7CB5Tk/inJRu2///6DXwjAAFXbqwMLVYW7qRpex52xAwCIhLADAIiEsAMAiISw\nAwCIhLADAIiEsAMAiISwAwCIhLADAIiEsAMAiISwAwCIhLADAIiEsAMAiISwAwCIhLADAIiE\nsAMAiISwAwCIhLADAIiEsAMAiISwAwCIhLADAIiEsAMAiISwAwCIhLADAIiEsAMAiISwAwCI\nhLADAIiEsAMAiISwAwCIhLADAIiEsAMAiISwAwCIhLADAIiEsAMAiISwAwCIhLADAIiEsAMA\niISwAwCIhLADAIiEsAMAiISwAwCIhLADAIiEsAMAiISwAwCIhLADAIiEsAMAiISwAwCIhLAD\nAIiEsAMAiISwAwCIhLADAIiEsAMAiISwAwCIhLADAIiEsAMAiISwAwCIhLADAIiEsAMAiESm\nXCsqPHDTJdff+YeZhZYtP3fIiQdum5RrxQAANaJMZ+zefuyyy3/xxkGnnHPulIPf+OV/Xfan\nt8uzXgCA2lGmsPv9tD9O+NqUiRuvt97GE6dM2ezRax8sz3oBAGpHmd6K/c309mPGDy9eHrH2\nbh3Trwxht75rzzvvvHvvvbd4uamp6Y477uh/aUmStLS0LKWhVtDo0aP7uTZJkhBCXV1d/7NV\nm6E1WmBJeY6XQWR3cpIkzc3Nzc3NlR5I6fWzp4ov4plMZvB7M5/P93NtmcLupc78GvVz15Wp\nXyPf+fK817a3t8+aNat4ube3t7jx/RvIPEPOAQccUOkhlF6UGwX08RwvA3fyUFGegEml+nu7\ntUxhl4RQ+MCED/y0/fbbr7XWWnMHlMm0tbX1v7TGxsaurq6enp5Bjuqaa64Z5BLKKZvN1tXV\n5fP5jo6OSo+lrJIkaWxsnDNnTqFQWPzcEWloaEilUt3d3V1dXZUeS1llMpm6urr29vZKD6Ss\nio/zEEJHR0f//x2PTy6XCyF0dnZWeiBllU6n6+vrQwi1eXDr7u4e/It4FeonYIov4r29vYM/\nuBUKhX7Od5Yp7FbOpV/s6BnfkAkh5Dv/nc6tPO+122yzzTbbbFO8XCgUpk+f3v/SGhoaOjs7\na+3VLp1Ol+oxMbSkUqnGxsb29vZaO/blcrlUKtXT01NrezyXy6XT6Vrb6r6w6+rqqsGDWwih\n1vZ4Npsthl1tHty6urpqLeVTqVQJX8T7Cbsy/fLEp8Y03PnP2cXLs/51Z8OY7cqzXgCA2lGm\nsNvyiM0fmnrZ4y+8/Mrzj1/2rd9POHLL8qwXAKB2lOmt2DEbH3fMdhdeeupxraFly88de+xG\nY8qzXgCA2lG2vzyRTDxkysRDyrU2AIDa42/FAgBEQtgBAERC2AEARELYAQBEQtgBAERC2AEA\nRELYAQBEQtgBAERC2AEARELYAQBEQtgBAERC2AEARELYAQBEQtgBAERC2AEARELYAQBEQtgB\nAERC2AEARELYAQBEQtgBAERC2AEARELYAQBEQtgBAERC2AEARELYAQBEQtgBAERC2AEARGJI\nht27777b1dVV6VGU26WXXjpx4sTJkydXeiDl1tvbO3369EKhUOmBlNvhhx8+ceLEq6++utID\nKbfOzs6ZM2dWehTl1tbWNnHixIkTJz700EOVHku5tbW1tbW1VXoU5fbQQw8V93gNbvvMmTM7\nOzsrPYpyu/rqqydOnHj44Ycv7RVllvYKllSSJGPGjKn0KKpRkiSzZ8/u7u52/9SIzs7O2bNn\np9Npe7wWtLa2zp49O4TQ1NRkj9eCpqam4h4fNWpUc3NzpYfDUpdOp2fPnt3Z2bm0n+BD8owd\nAAALEnYAAJGourdiWZR11lln0qRJq622WqUHQplMmDBhtdVWGzduXKUHQjlkMplJkyaFELwP\nWyPGjBlT3OOZjBfimjBu3LhJkyYtt9xyS3tFSQ1+Jh0AIEreigUAiISwAwCIhLADAIiEz2xW\n3ltP3HXlD+546oVXe7LDPrblTl/58r4jM0kI4eID93pg5vtf4bjshPOuOW39EEIIhQduuuT6\nO/8ws9Cy5ecOOfHAbZO5syxqOlWkdLvV7h4adtlllwUn3nHHHZ7gUbr0oL2/euOP5plQque1\nvV+l5tvj1fCCLuwqrPPd3x171ndX3Hb/M476eF3763fdcNWUS1a9+pQtQ6H7j61dmx1z4qQR\nueKc9aNWLV54+7HLLv/FG8eecs6y4Y3vf/u/Llv3o8dvOqaf6VSR0u1Wu3uoOO200+b98fk7\nrro3/VlP8NgUut954+U/3fs/983o+Oo8k0v1vLb3q87C9ni1vKAXqKjnrj9uz32/3tM798fu\nOc/suschhUKhc9ajO++881Nt3Qve5GdH73/hX94uXn7rz1P3O/qn/U+nepRwt9rdQ1Hry3fv\ns8+JL3f2eIJH5t1nLj388MMPP/zwnXfeed7ppXpe2/vVZqF7vEpe0H3GrsKaV9luv/0PSr93\njjVJUklSF0LonPFQkqq//+Kv7bvnbp/f77ALpv2i470vpvnN9PZdxg8vXh6x9m4d0+/rfzrV\no4S71e4ecgq9cy477drtv37aitm0J3hkRqx9/LRp06ZNmzbf9FI9r+39arPQPV4lL+jeiq2w\n5Sbuuud7l/Ptr11/7rmr73RSCGHGky+HJF2/8e7nH7xy+xvP3XrV905Nr3HpoeuGEF7qzK9R\nP3fHZerXyHe+XLy8qOlUjxLuVrt7yHnhJ998evTnT1t/VPAErxmlel7b+0NClbygC7uq0Nv9\nzq9vu+XWn/523A6HX3jo+iGEFSaeddPWdS3N2RBCWGWVU1ftPOC4aeHQi0MISQgf/FLpuT8t\najrVo4S71e4eWvKdL533w2f3v+obxR89wWtEqZ7X9v4QUvEXdG/FVt70J+864dCjfvNq04mX\nf//Uw7bPJCGEkK5vmvsgCCGE0DB2Ur7jxeLllXPpFzt6ipfznf9O51bufzrVo4S71e4eWl68\n/ZKO0XvsuExD8UdP8BpRque1vT9UVMMLurCrsI637jvuzOu2mHzZ1BMPWW+5xr7pPzntq6df\n9Zf3Z5v+QLZl0+LlT41puPOfs4uXZ/3rzoYx2/U/nepRwt1qdw8lhfzV//ufdQ7boW+CJ3iN\nKNXz2t4fEqrkBT39jW98YxBbwWD9aeo5v+/ccsdNml6ex0orrTSi8MR1N98yq2HZYZn8K889\nfvXUa1ba/5Stxw0LIYxe/pXvXnT/uA3G9b7z/JXn3Tz+y1/afPnGfqZTPUq4W+3uIWTOW7d/\n7/a/nHLcoSMzc/8v7Qkeq1tuuWW//fbr+7FUz2t7v2rNu8er5AU9KRS8VV9J831pYdEdd9wR\nQvjLL2+46ee/ef7VGbkRK2y7yyFH77l53/cW3nf9hdf98uHW0LLl5w4/6aBPLm46VaR0u9Xu\nHjJe+MHxU+4ec9v/nD7vRE/wKO2yyy7FA/h7SvW8tver1Lx7vEpe0IUdAEAkfMYOACASwg4A\nIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOACASwg4AIBLCDgAgEsIOqEaPnbph0q+Nv/nn0q6x\nuNj5fpxXOts0fpNJF9xa4vX2PwyAJZKp9AAAFmLZrQ86sfvVvh8vvvjiEMKJJ57YN2XFLZZd\n7EJ22GGHEMLdd9/9oYcxzxoLs1575qZb7/ravhu/Ouyly3Zc6UMvE2Dp8bdigSGgeBJrSY9X\nS3Sr+WZe6G2fuWHndQ+5s2XFybNevmyJRjJwH25LAYqEHTAEVEnY9bQ/V9c4PpUZnu+esUQj\nGThhBwyGz9gBQ1t369PnfWmvdVYck0tnRi4/bo+jznpqVld4r5DCvJ9a6+247TtTJm6yzoim\nXCpdN3qFcbscdtrfZnUNfF2Z+lWLi/zA1ELn7ZeesOW6qzRl08OWWX33L37z2dbuuVf1u8ZC\nzzvfP/eYCeuu3JzLNI5Yfru9v/Lwm+3zLDd/51Wnf/rj6w1vqk9n6lcYt+Fhp175bo/gA/rj\njB0wBCzqPFa+8z+fHfeRX73cut5Oh+y31bh//eHW6+78W+Nyn3riX//3/P337LjjjiGEu+66\nK4Swww47/OSg9T7/g6fT2WUOPu6L40clD99+5c8ee3v0R059+8nzF1zFQtfYMf1nDWN2HzHu\n9HefO6dv4qV7r33Cj/+xylafP3THjWc8e/8VP7infqUd//7cHSvn0v2ssZCfefSEta7+01vN\nq044dL/PjEnPvu/H1z705jqdMx4prvee4zfb/rI/tqy29ZFf+PTYXOcDP7zy7qfeXX3P//nn\nbfsvtbsZGPoKAFVvUcerR762QQhhtd2/1zflur3WCCGsf8JDC95q2Ww6hDDhv58q/tjT/lwI\nIUnVLXQVC65x9uvPnLHjyuncitc+/W7fxNd+96UQwkqfuai7d+6UJ288KISw0Wl/7H+Nf/32\nliGEkese9UZXvnhtb8+M4zcc07feVXKZEMJjs7uK13a1/nXZZZddYeX1l/CeA2qLM3bAELCo\nM3b7L9t0y5tzrnil9ZgVmopT2l77bvMKX2ocu3fbm7fOd6u33347hJAbObolnYQQOqY/1DBm\ny74ZFnrGbv5hpOqOu/Gxyw5Yv2/KVR8be8zf3r7ildnHrND83rTez41uvj+7X+tr1/azxkOW\na77hjbaznn33G+NH9C1t5gvnjBh3ZnGG9ZqyT8/p/uqND3zrC9vkfP8JMDDCDhgCFhV2o7Pp\nd7p73+3uHZGZ2z6F/MxUZkQqMyzfPXPBW73xxD033n7Pk8/844UXXnjyib/P7OkN/Ybd008/\n3Xfbzlmv3Xj6gZfc88puFzz201M2Lk7cqCX3l9auR5/8ezHdiu7+3BYnvdzU0/FyP2scVZd+\nt6d3enfvqMz7N+ztmZ6uG1Oc4dmbTt76sEvf6s5nh68w4RNbb7X1Np/bc/8t1nq/AgEWJOyA\nIWBRYdeSSbfm5wu72anMsFS6Kd/TOt+t/nblQRsde1OqcbW9D9rr4xt+dKPNJn5yo5VDv2E3\n3xo73r27YdSOuWGf6Jj5YHFKQzrV0buQo2iSpHt7e/pZ48i69Iye3je78mPr3v8ltkLPu6m6\nUX3r7Z71r5//5Of3//a3v/3tb5/451tJKvuFix+68fiNB3FHApETdsAQsKiw23Ns0+1vz7nq\n1dajl3/vrdjXpzUvf2TDmN3nvHX7fLcqnl279aXZe69UfNu0N0nSYUnCLoR8kmSSpK63d+5v\ntq7TlH12Tvd8fdannzUesXzzta+3ff3pd85dZ2Tf/DP/ed6INU8vzvDggw9m6lefsOmKxaue\nv/fCtSZNydSv0d3+wpLee0Dt8HUnwBB20hfWCCFcOPmHfVNunXx+CGHNA05ecOZ/dvSEEDYb\n21D88dX7z1ryFRazrLvvW0fO3GLZEMJhVz/eN0fnu4+uv+yYdT/7k/7X+JWTNgwhXLn7lDe6\ne4tTCvlZ39zrO30zHLb9p7bd7vMz3lvTqpt/OoQQQu+SjxmoIc7YAUPAos7Y9XQ8/5nVNvjN\nG3PW/9zh+2y15ou//+G0n/+1cblJT/zr7nH16eL7pF+56KoNR2UOOfSI/952xeMeeHW1Hb40\n5fMbT3/6gUu/94dc979f7cyfMe2WQw/Ya42GTFj8GbtQfAv16Tnd6zRkQggd0+/9+Go7PdXW\nPWG3wz+7+fiut5//8bQb/jFnxE3PP7/vys39rPGQfbc/f+u1r/3z282rTDhkv+2Xzc2570fT\nHi/sMuPZG4vrveuIj+x07d9HfmTS4XtOHJtr//UN/3XPP2ZuPPn/Hrvs00v5zgaGsqX8W7cA\nJdDP8apzxt/OOGKX1cYOS5Jk2DJr7HbkmU/NnPsVIb86++BVxzSFEBqGjSkUCt1zXjjtoM8s\nP6IhSeXW3WqvO56Z8fvz9x/bmElnh//orTnzrWJRa9x/mcYQwmd/8FzflDmvP/K1w3ZZc7mR\nqSSpH7bMVrscfseT7xSv6n+N+e63vnfGURuPWz6bSuqax263/5Rn2rr71tubb73mjCM2GrdC\nLpNKpetXWW/CVy+6rat3wREBvM8ZOwCASPiMHQBAJIQdAEAkhB0AQCSEHQBAJIQdAEAkhB0A\nQCSEHQBAJIQdAEAkhB0AQCSEHQBAJIQdAEAkhB0AQCSEHQBAJP4f1gqhxpD1xY4AAAAASUVO\nRK5CYII=",
      "text/plain": [
       "plot without title"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "alz_sdt = data.table(as(sample_data(alz_nt), \"data.frame\"),\n",
    "                 TotalReads = sample_sums(alz_nt), keep.rownames = TRUE)\n",
    "setnames(alz_sdt, \"rn\", \"SampleName\")\n",
    "pSeqDepth = ggplot(alz_sdt, aes(TotalReads)) + geom_histogram() + ggtitle(\"Sequencing Depth\")\n",
    "pSeqDepth\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": true
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "You set `rngseed` to FALSE. Make sure you've set & recorded\n",
      " the random seed of your session for reproducibility.\n",
      "See `?set.seed`\n",
      "\n",
      "...\n",
      "622OTUs were removed because they are no longer \n",
      "present in any sample after random subsampling\n",
      "\n",
      "...\n",
      "You set `rngseed` to FALSE. Make sure you've set & recorded\n",
      " the random seed of your session for reproducibility.\n",
      "See `?set.seed`\n",
      "\n",
      "...\n",
      "9 samples removedbecause they contained fewer reads than `sample.size`.\n",
      "Up to first five removed samples are: \n",
      "\n",
      "M10M36M11M34M23\t\n",
      "...\n",
      "217OTUs were removed because they are no longer \n",
      "present in any sample after random subsampling\n",
      "\n",
      "...\n",
      "You set `rngseed` to FALSE. Make sure you've set & recorded\n",
      " the random seed of your session for reproducibility.\n",
      "See `?set.seed`\n",
      "\n",
      "...\n",
      "15 samples removedbecause they contained fewer reads than `sample.size`.\n",
      "Up to first five removed samples are: \n",
      "\n",
      "M29M26M10M21M36\t\n",
      "...\n",
      "415OTUs were removed because they are no longer \n",
      "present in any sample after random subsampling\n",
      "\n",
      "...\n"
     ]
    }
   ],
   "source": [
    "#Rarefy all data to 1k\n",
    "alz_nt_1k = rarefy_even_depth(alz_nt, sample.size = 1000,\n",
    "                                  replace = FALSE, trimOTUs = TRUE, verbose = TRUE)\n",
    "alz_t_1k <- alz_nt_1k %>%\n",
    "    transform_sample_counts(function(x) x / sum(x))\n",
    "saveRDS(alz_nt_1k, \"/gpfs/data/gilbert-lab/Lutz/Alzheimers/RDSfiles/alz_nt_1k.rds\")                                                   \n",
    "saveRDS(alz_t_1k, \"/gpfs/data/gilbert-lab/Lutz/Alzheimers/RDSfiles/alz_t_1k.rds\")\n",
    "\n",
    "#Rarefy all data to 5k\n",
    "alz_nt_5k = rarefy_even_depth(alz_nt, sample.size = 5000,\n",
    "                                 replace = FALSE, trimOTUs = TRUE, verbose = TRUE)\n",
    "alz_t_5k <- alz_nt_5k %>%\n",
    "    transform_sample_counts(function(x) x / sum(x))\n",
    "saveRDS(alz_nt_5k, \"/gpfs/data/gilbert-lab/Lutz/Alzheimers/RDSfiles/alz_nt_5k.rds\")                                                   \n",
    "saveRDS(alz_t_5k, \"/gpfs/data/gilbert-lab/Lutz/Alzheimers/RDSfiles/alz_t_5k.rds\")\n",
    "                            \n",
    "#Rarefy all data to 10k\n",
    "alz_nt_10k = rarefy_even_depth(alz_nt, sample.size = 10000,\n",
    "                                  replace = FALSE, trimOTUs = TRUE, verbose = TRUE)\n",
    "alz_t_10k <- alz_nt_10k %>%\n",
    "    transform_sample_counts(function(x) x / sum(x))\n",
    "saveRDS(alz_nt_10k, \"/gpfs/data/gilbert-lab/Lutz/Alzheimers/RDSfiles/alz_nt_10k.rds\")                                                   \n",
    "saveRDS(alz_t_10k, \"/gpfs/data/gilbert-lab/Lutz/Alzheimers/RDSfiles/alz_t_10k.rds\")\n",
    "                            "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "phyloseq-class experiment-level object\n",
       "otu_table()   OTU Table:         [ 3441 taxa and 49 samples ]\n",
       "sample_data() Sample Data:       [ 49 samples by 12 sample variables ]\n",
       "tax_table()   Taxonomy Table:    [ 3441 taxa by 7 taxonomic ranks ]\n",
       "phy_tree()    Phylogenetic Tree: [ 3441 tips and 3439 internal nodes ]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "phyloseq-class experiment-level object\n",
       "otu_table()   OTU Table:         [ 3430 taxa and 49 samples ]\n",
       "sample_data() Sample Data:       [ 49 samples by 12 sample variables ]\n",
       "tax_table()   Taxonomy Table:    [ 3430 taxa by 7 taxonomic ranks ]\n",
       "phy_tree()    Phylogenetic Tree: [ 3430 tips and 3428 internal nodes ]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "phyloseq-class experiment-level object\n",
       "otu_table()   OTU Table:         [ 2819 taxa and 49 samples ]\n",
       "sample_data() Sample Data:       [ 49 samples by 12 sample variables ]\n",
       "tax_table()   Taxonomy Table:    [ 2819 taxa by 7 taxonomic ranks ]\n",
       "phy_tree()    Phylogenetic Tree: [ 2819 tips and 2817 internal nodes ]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "phyloseq-class experiment-level object\n",
       "otu_table()   OTU Table:         [ 3224 taxa and 40 samples ]\n",
       "sample_data() Sample Data:       [ 40 samples by 12 sample variables ]\n",
       "tax_table()   Taxonomy Table:    [ 3224 taxa by 7 taxonomic ranks ]\n",
       "phy_tree()    Phylogenetic Tree: [ 3224 tips and 3222 internal nodes ]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "phyloseq-class experiment-level object\n",
       "otu_table()   OTU Table:         [ 3026 taxa and 34 samples ]\n",
       "sample_data() Sample Data:       [ 34 samples by 12 sample variables ]\n",
       "tax_table()   Taxonomy Table:    [ 3026 taxa by 7 taxonomic ranks ]\n",
       "phy_tree()    Phylogenetic Tree: [ 3026 tips and 3024 internal nodes ]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "alz_nt\n",
    "alz\n",
    "alz_nt_1k\n",
    "alz_nt_5k\n",
    "alz_nt_10k"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 86,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th></th><th scope=col>SampleID</th><th scope=col>BarcodeSequence</th><th scope=col>LinkerPrimerSequence</th><th scope=col>PrimerPlate</th><th scope=col>PrimerWell</th><th scope=col>Tube</th><th scope=col>MouseTag</th><th scope=col>SampleType</th><th scope=col>Sex</th><th scope=col>Treatment</th><th scope=col>FecalMassGm</th><th scope=col>Description</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><th scope=row>X17</th><td>X17                  </td><td>GCATATGCACTG         </td><td>GTGTGYCAGCMGCCGCGGTAA</td><td>D5                   </td><td>ABX_cecal_plate1     </td><td>3046                 </td><td>M                    </td><td>Cecal                </td><td>Abx_FMT              </td><td>M_Abx_FMT            </td><td>0.0493               </td><td>M17                  </td></tr>\n",
       "\t<tr><th scope=row>X5.F.7.Abx</th><td>X5.F.7.Abx           </td><td>AATAACTAGGGT         </td><td>GTGTGYCAGCMGCCGCGGTAA</td><td>D11                  </td><td>ABX_fecal_plate2     </td><td>unknown              </td><td>F                    </td><td>Fecal                </td><td>Abx                  </td><td>F_Abx                </td><td>Unknown              </td><td>5.F.7.Abx            </td></tr>\n",
       "\t<tr><th scope=row>X2.M.1.Abx</th><td>X2.M.1.Abx           </td><td>ACTACTGAGGAT         </td><td>GTGTGYCAGCMGCCGCGGTAA</td><td>A11                  </td><td>ABX_fecal_plate2     </td><td>unknown              </td><td>M                    </td><td>Fecal                </td><td>Abx                  </td><td>M_Abx                </td><td>Unknown              </td><td>2.M.1.Abx            </td></tr>\n",
       "\t<tr><th scope=row>X9</th><td>X9                   </td><td>CATTTCGCACTT         </td><td>GTGTGYCAGCMGCCGCGGTAA</td><td>C9                   </td><td>ABX_cecal_plate1     </td><td>3190                 </td><td>M                    </td><td>Cecal                </td><td>Control              </td><td>M_Control            </td><td>0.0299               </td><td>M9                   </td></tr>\n",
       "\t<tr><th scope=row>X4.F.7.ctr</th><td>X4.F.7.ctr           </td><td>ACCTTACACCTT         </td><td>GTGTGYCAGCMGCCGCGGTAA</td><td>D2                   </td><td>ABX_fecal_plate2     </td><td>unknown              </td><td>F                    </td><td>Fecal                </td><td>Control              </td><td>F_Control            </td><td>Unknown              </td><td>4.F.7.ctr            </td></tr>\n",
       "\t<tr><th scope=row>X4.F.1.Ctr</th><td>X4.F.1.Ctr           </td><td>AAGAACTCATGA         </td><td>GTGTGYCAGCMGCCGCGGTAA</td><td>C6                   </td><td>ABX_fecal_plate2     </td><td>unknown              </td><td>F                    </td><td>Fecal                </td><td>Control              </td><td>F_Control            </td><td>Unknown              </td><td>4.F.1.Ctr            </td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|llllllllllll}\n",
       "  & SampleID & BarcodeSequence & LinkerPrimerSequence & PrimerPlate & PrimerWell & Tube & MouseTag & SampleType & Sex & Treatment & FecalMassGm & Description\\\\\n",
       "\\hline\n",
       "\tX17 & X17                   & GCATATGCACTG          & GTGTGYCAGCMGCCGCGGTAA & D5                    & ABX\\_cecal\\_plate1  & 3046                  & M                     & Cecal                 & Abx\\_FMT             & M\\_Abx\\_FMT         & 0.0493                & M17                  \\\\\n",
       "\tX5.F.7.Abx & X5.F.7.Abx            & AATAACTAGGGT          & GTGTGYCAGCMGCCGCGGTAA & D11                   & ABX\\_fecal\\_plate2  & unknown               & F                     & Fecal                 & Abx                   & F\\_Abx               & Unknown               & 5.F.7.Abx            \\\\\n",
       "\tX2.M.1.Abx & X2.M.1.Abx            & ACTACTGAGGAT          & GTGTGYCAGCMGCCGCGGTAA & A11                   & ABX\\_fecal\\_plate2  & unknown               & M                     & Fecal                 & Abx                   & M\\_Abx               & Unknown               & 2.M.1.Abx            \\\\\n",
       "\tX9 & X9                    & CATTTCGCACTT          & GTGTGYCAGCMGCCGCGGTAA & C9                    & ABX\\_cecal\\_plate1  & 3190                  & M                     & Cecal                 & Control               & M\\_Control           & 0.0299                & M9                   \\\\\n",
       "\tX4.F.7.ctr & X4.F.7.ctr            & ACCTTACACCTT          & GTGTGYCAGCMGCCGCGGTAA & D2                    & ABX\\_fecal\\_plate2  & unknown               & F                     & Fecal                 & Control               & F\\_Control           & Unknown               & 4.F.7.ctr            \\\\\n",
       "\tX4.F.1.Ctr & X4.F.1.Ctr            & AAGAACTCATGA          & GTGTGYCAGCMGCCGCGGTAA & C6                    & ABX\\_fecal\\_plate2  & unknown               & F                     & Fecal                 & Control               & F\\_Control           & Unknown               & 4.F.1.Ctr            \\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "| <!--/--> | SampleID | BarcodeSequence | LinkerPrimerSequence | PrimerPlate | PrimerWell | Tube | MouseTag | SampleType | Sex | Treatment | FecalMassGm | Description |\n",
       "|---|---|---|---|---|---|---|---|---|---|---|---|---|\n",
       "| X17 | X17                   | GCATATGCACTG          | GTGTGYCAGCMGCCGCGGTAA | D5                    | ABX_cecal_plate1      | 3046                  | M                     | Cecal                 | Abx_FMT               | M_Abx_FMT             | 0.0493                | M17                   |\n",
       "| X5.F.7.Abx | X5.F.7.Abx            | AATAACTAGGGT          | GTGTGYCAGCMGCCGCGGTAA | D11                   | ABX_fecal_plate2      | unknown               | F                     | Fecal                 | Abx                   | F_Abx                 | Unknown               | 5.F.7.Abx             |\n",
       "| X2.M.1.Abx | X2.M.1.Abx            | ACTACTGAGGAT          | GTGTGYCAGCMGCCGCGGTAA | A11                   | ABX_fecal_plate2      | unknown               | M                     | Fecal                 | Abx                   | M_Abx                 | Unknown               | 2.M.1.Abx             |\n",
       "| X9 | X9                    | CATTTCGCACTT          | GTGTGYCAGCMGCCGCGGTAA | C9                    | ABX_cecal_plate1      | 3190                  | M                     | Cecal                 | Control               | M_Control             | 0.0299                | M9                    |\n",
       "| X4.F.7.ctr | X4.F.7.ctr            | ACCTTACACCTT          | GTGTGYCAGCMGCCGCGGTAA | D2                    | ABX_fecal_plate2      | unknown               | F                     | Fecal                 | Control               | F_Control             | Unknown               | 4.F.7.ctr             |\n",
       "| X4.F.1.Ctr | X4.F.1.Ctr            | AAGAACTCATGA          | GTGTGYCAGCMGCCGCGGTAA | C6                    | ABX_fecal_plate2      | unknown               | F                     | Fecal                 | Control               | F_Control             | Unknown               | 4.F.1.Ctr             |\n",
       "\n"
      ],
      "text/plain": [
       "           SampleID   BarcodeSequence LinkerPrimerSequence  PrimerPlate\n",
       "X17        X17        GCATATGCACTG    GTGTGYCAGCMGCCGCGGTAA D5         \n",
       "X5.F.7.Abx X5.F.7.Abx AATAACTAGGGT    GTGTGYCAGCMGCCGCGGTAA D11        \n",
       "X2.M.1.Abx X2.M.1.Abx ACTACTGAGGAT    GTGTGYCAGCMGCCGCGGTAA A11        \n",
       "X9         X9         CATTTCGCACTT    GTGTGYCAGCMGCCGCGGTAA C9         \n",
       "X4.F.7.ctr X4.F.7.ctr ACCTTACACCTT    GTGTGYCAGCMGCCGCGGTAA D2         \n",
       "X4.F.1.Ctr X4.F.1.Ctr AAGAACTCATGA    GTGTGYCAGCMGCCGCGGTAA C6         \n",
       "           PrimerWell       Tube    MouseTag SampleType Sex     Treatment\n",
       "X17        ABX_cecal_plate1 3046    M        Cecal      Abx_FMT M_Abx_FMT\n",
       "X5.F.7.Abx ABX_fecal_plate2 unknown F        Fecal      Abx     F_Abx    \n",
       "X2.M.1.Abx ABX_fecal_plate2 unknown M        Fecal      Abx     M_Abx    \n",
       "X9         ABX_cecal_plate1 3190    M        Cecal      Control M_Control\n",
       "X4.F.7.ctr ABX_fecal_plate2 unknown F        Fecal      Control F_Control\n",
       "X4.F.1.Ctr ABX_fecal_plate2 unknown F        Fecal      Control F_Control\n",
       "           FecalMassGm Description\n",
       "X17        0.0493      M17        \n",
       "X5.F.7.Abx Unknown     5.F.7.Abx  \n",
       "X2.M.1.Abx Unknown     2.M.1.Abx  \n",
       "X9         0.0299      M9         \n",
       "X4.F.7.ctr Unknown     4.F.7.ctr  \n",
       "X4.F.1.Ctr Unknown     4.F.1.Ctr  "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "sample_tab = data.frame(sample_data(alz_nt))\n",
    "names(sample_tab) <- c(\"SampleID\", \"BarcodeSequence\", \"LinkerPrimerSequence\", \"PrimerPlate\", \"PrimerWell\", \"Tube\", \"MouseTag\", \"SampleType\", \"Sex\", \"Treatment\", \"FecalMassGm\",\"Description\")\n",
    "head(sample_tab)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 87,
   "metadata": {},
   "outputs": [],
   "source": [
    "alz_nt_fecal = subset_samples(alz_nt, SampleType==\"Fecal\")\n",
    "alz_nt_cecal = subset_samples(alz_nt, SampleType==\"Cecal\")\n",
    "\n",
    "alz_nt_fec_M = subset_samples(alz_nt_fecal, Sex==\"M\")\n",
    "alz_nt_fec_F = subset_samples(alz_nt_fecal, Sex==\"F\")\n",
    "alz_nt_cec_M = subset_samples(alz_nt_cecal, Sex==\"M\")\n",
    "alz_nt_cec_F = subset_samples(alz_nt_cecal, Sex==\"F\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 88,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th scope=col>Treatment</th><th scope=col>median</th><th scope=col>mean</th><th scope=col>sd</th><th scope=col>sem</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><td>M_Abx    </td><td>163      </td><td>166.6000 </td><td>21.44865 </td><td> 6.782658</td></tr>\n",
       "\t<tr><td>M_Abx_FMT</td><td>185      </td><td>181.1111 </td><td>26.01175 </td><td> 8.670583</td></tr>\n",
       "\t<tr><td>M_Control</td><td>203      </td><td>196.9000 </td><td>47.94082 </td><td>15.160218</td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|lllll}\n",
       " Treatment & median & mean & sd & sem\\\\\n",
       "\\hline\n",
       "\t M\\_Abx     & 163         & 166.6000    & 21.44865    &  6.782658  \\\\\n",
       "\t M\\_Abx\\_FMT & 185           & 181.1111      & 26.01175      &  8.670583    \\\\\n",
       "\t M\\_Control & 203         & 196.9000    & 47.94082    & 15.160218  \\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "| Treatment | median | mean | sd | sem |\n",
       "|---|---|---|---|---|\n",
       "| M_Abx     | 163       | 166.6000  | 21.44865  |  6.782658 |\n",
       "| M_Abx_FMT | 185       | 181.1111  | 26.01175  |  8.670583 |\n",
       "| M_Control | 203       | 196.9000  | 47.94082  | 15.160218 |\n",
       "\n"
      ],
      "text/plain": [
       "  Treatment median mean     sd       sem      \n",
       "1 M_Abx     163    166.6000 21.44865  6.782658\n",
       "2 M_Abx_FMT 185    181.1111 26.01175  8.670583\n",
       "3 M_Control 203    196.9000 47.94082 15.160218"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "\n",
       "\tKruskal-Wallis rank sum test\n",
       "\n",
       "data:  Observed by Treatment\n",
       "Kruskal-Wallis chi-squared = 6.7183, df = 2, p-value = 0.03476\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  Kruskal-Wallis rank sum test\n",
      "\n",
      "data: x and group\n",
      "Kruskal-Wallis chi-squared = 6.7183, df = 2, p-value = 0.03\n",
      "\n",
      "\n",
      "                           Comparison of x by group                            \n",
      "                                 (Bonferroni)                                  \n",
      "Col Mean-|\n",
      "Row Mean |      M_Abx   M_Abx_FM\n",
      "---------+----------------------\n",
      "M_Abx_FM |  -1.363169\n",
      "         |     0.2592\n",
      "         |\n",
      "M_Contro |  -2.589289  -1.157059\n",
      "         |    0.0144*     0.3709\n",
      "\n",
      "alpha = 0.05\n",
      "Reject Ho if p <= alpha/2\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<dl>\n",
       "\t<dt>$chi2</dt>\n",
       "\t\t<dd>6.71833388175936</dd>\n",
       "\t<dt>$Z</dt>\n",
       "\t\t<dd><ol class=list-inline>\n",
       "\t<li>-1.36316947170347</li>\n",
       "\t<li>-2.58928900889213</li>\n",
       "\t<li>-1.157059384741</li>\n",
       "</ol>\n",
       "</dd>\n",
       "\t<dt>$P</dt>\n",
       "\t\t<dd><ol class=list-inline>\n",
       "\t<li>0.0864145530663608</li>\n",
       "\t<li>0.00480871698658505</li>\n",
       "\t<li>0.123624049468183</li>\n",
       "</ol>\n",
       "</dd>\n",
       "\t<dt>$P.adjusted</dt>\n",
       "\t\t<dd><ol class=list-inline>\n",
       "\t<li>0.259243659199082</li>\n",
       "\t<li>0.0144261509597551</li>\n",
       "\t<li>0.37087214840455</li>\n",
       "</ol>\n",
       "</dd>\n",
       "\t<dt>$comparisons</dt>\n",
       "\t\t<dd><ol class=list-inline>\n",
       "\t<li>'M_Abx - M_Abx_FMT'</li>\n",
       "\t<li>'M_Abx - M_Control'</li>\n",
       "\t<li>'M_Abx_FMT - M_Control'</li>\n",
       "</ol>\n",
       "</dd>\n",
       "</dl>\n"
      ],
      "text/latex": [
       "\\begin{description}\n",
       "\\item[\\$chi2] 6.71833388175936\n",
       "\\item[\\$Z] \\begin{enumerate*}\n",
       "\\item -1.36316947170347\n",
       "\\item -2.58928900889213\n",
       "\\item -1.157059384741\n",
       "\\end{enumerate*}\n",
       "\n",
       "\\item[\\$P] \\begin{enumerate*}\n",
       "\\item 0.0864145530663608\n",
       "\\item 0.00480871698658505\n",
       "\\item 0.123624049468183\n",
       "\\end{enumerate*}\n",
       "\n",
       "\\item[\\$P.adjusted] \\begin{enumerate*}\n",
       "\\item 0.259243659199082\n",
       "\\item 0.0144261509597551\n",
       "\\item 0.37087214840455\n",
       "\\end{enumerate*}\n",
       "\n",
       "\\item[\\$comparisons] \\begin{enumerate*}\n",
       "\\item 'M\\_Abx - M\\_Abx\\_FMT'\n",
       "\\item 'M\\_Abx - M\\_Control'\n",
       "\\item 'M\\_Abx\\_FMT - M\\_Control'\n",
       "\\end{enumerate*}\n",
       "\n",
       "\\end{description}\n"
      ],
      "text/markdown": [
       "$chi2\n",
       ":   6.71833388175936\n",
       "$Z\n",
       ":   1. -1.36316947170347\n",
       "2. -2.58928900889213\n",
       "3. -1.157059384741\n",
       "\n",
       "\n",
       "\n",
       "$P\n",
       ":   1. 0.0864145530663608\n",
       "2. 0.00480871698658505\n",
       "3. 0.123624049468183\n",
       "\n",
       "\n",
       "\n",
       "$P.adjusted\n",
       ":   1. 0.259243659199082\n",
       "2. 0.0144261509597551\n",
       "3. 0.37087214840455\n",
       "\n",
       "\n",
       "\n",
       "$comparisons\n",
       ":   1. 'M_Abx - M_Abx_FMT'\n",
       "2. 'M_Abx - M_Control'\n",
       "3. 'M_Abx_FMT - M_Control'\n",
       "\n",
       "\n",
       "\n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "$chi2\n",
       "[1] 6.718334\n",
       "\n",
       "$Z\n",
       "[1] -1.363169 -2.589289 -1.157059\n",
       "\n",
       "$P\n",
       "[1] 0.086414553 0.004808717 0.123624049\n",
       "\n",
       "$P.adjusted\n",
       "[1] 0.25924366 0.01442615 0.37087215\n",
       "\n",
       "$comparisons\n",
       "[1] \"M_Abx - M_Abx_FMT\"     \"M_Abx - M_Control\"     \"M_Abx_FMT - M_Control\"\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "#Observed richness by treatment group - MALES - FECAL\n",
    "\n",
    "obs = estimate_richness(na.omit(alz_nt_fec_M), split = TRUE, measures = \"Observed\")\n",
    "obs$SampleID = rownames(obs)\n",
    "\n",
    "dfobs = merge(obs, sample_tab, on=\"SampleID\")\n",
    "\n",
    "ordobs = obs %>% group_by(\"Treatment\") %>%\n",
    "    summarize(med.obs = median(Observed)) %>%\n",
    "    arrange(desc(med.obs)) %>%\n",
    "    magrittr::extract2(\"Treatment\")\n",
    "\n",
    "stats.obs = dfobs %>%\n",
    "    group_by(Treatment) %>%\n",
    "    summarize(median = median(Observed), \n",
    "              mean = mean(Observed),\n",
    "              sd = sd(Observed),\n",
    "              sem = sd(Observed/sqrt(length(Observed)))) #%>%\n",
    "    #arrange(-desc(mean))\n",
    "stats.obs\n",
    "\n",
    "#CHECK FOR SIGNIFICANCE\n",
    "stats.krusk = kruskal.test(Observed~Treatment, data=dfobs)\n",
    "stats.krusk\n",
    "\n",
    "#POST-HOC ANALYSIS\n",
    "dunn = dunn.test(dfobs$Observed, g=dfobs$Treatment, method=\"bonferroni\")\n",
    "dunn"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 89,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th scope=col>Treatment</th><th scope=col>median</th><th scope=col>mean</th><th scope=col>sd</th><th scope=col>sem</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><td>M_Abx    </td><td>139      </td><td>150.5000 </td><td>47.72665 </td><td>15.09249 </td></tr>\n",
       "\t<tr><td>M_Abx_FMT</td><td>158      </td><td>174.3333 </td><td>55.82114 </td><td>18.60705 </td></tr>\n",
       "\t<tr><td>M_Control</td><td>191      </td><td>191.8000 </td><td>55.81278 </td><td>17.64955 </td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|lllll}\n",
       " Treatment & median & mean & sd & sem\\\\\n",
       "\\hline\n",
       "\t M\\_Abx     & 139         & 150.5000    & 47.72665    & 15.09249   \\\\\n",
       "\t M\\_Abx\\_FMT & 158           & 174.3333      & 55.82114      & 18.60705     \\\\\n",
       "\t M\\_Control & 191         & 191.8000    & 55.81278    & 17.64955   \\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "| Treatment | median | mean | sd | sem |\n",
       "|---|---|---|---|---|\n",
       "| M_Abx     | 139       | 150.5000  | 47.72665  | 15.09249  |\n",
       "| M_Abx_FMT | 158       | 174.3333  | 55.82114  | 18.60705  |\n",
       "| M_Control | 191       | 191.8000  | 55.81278  | 17.64955  |\n",
       "\n"
      ],
      "text/plain": [
       "  Treatment median mean     sd       sem     \n",
       "1 M_Abx     139    150.5000 47.72665 15.09249\n",
       "2 M_Abx_FMT 158    174.3333 55.82114 18.60705\n",
       "3 M_Control 191    191.8000 55.81278 17.64955"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "\n",
       "\tKruskal-Wallis rank sum test\n",
       "\n",
       "data:  Observed by Treatment\n",
       "Kruskal-Wallis chi-squared = 2.9875, df = 2, p-value = 0.2245\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  Kruskal-Wallis rank sum test\n",
      "\n",
      "data: x and group\n",
      "Kruskal-Wallis chi-squared = 2.9875, df = 2, p-value = 0.22\n",
      "\n",
      "\n",
      "                           Comparison of x by group                            \n",
      "                                 (Bonferroni)                                  \n",
      "Col Mean-|\n",
      "Row Mean |      M_Abx   M_Abx_FM\n",
      "---------+----------------------\n",
      "M_Abx_FM |  -0.981494\n",
      "         |     0.4895\n",
      "         |\n",
      "M_Contro |  -1.720538  -0.693154\n",
      "         |     0.1280     0.7323\n",
      "\n",
      "alpha = 0.05\n",
      "Reject Ho if p <= alpha/2\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<dl>\n",
       "\t<dt>$chi2</dt>\n",
       "\t\t<dd>2.98748699414051</dd>\n",
       "\t<dt>$Z</dt>\n",
       "\t\t<dd><ol class=list-inline>\n",
       "\t<li>-0.981494767538991</li>\n",
       "\t<li>-1.72053810052189</li>\n",
       "\t<li>-0.693154047118709</li>\n",
       "</ol>\n",
       "</dd>\n",
       "\t<dt>$P</dt>\n",
       "\t\t<dd><ol class=list-inline>\n",
       "\t<li>0.163174406730601</li>\n",
       "\t<li>0.0426673370880744</li>\n",
       "\t<li>0.244106441421366</li>\n",
       "</ol>\n",
       "</dd>\n",
       "\t<dt>$P.adjusted</dt>\n",
       "\t\t<dd><ol class=list-inline>\n",
       "\t<li>0.489523220191803</li>\n",
       "\t<li>0.128002011264223</li>\n",
       "\t<li>0.732319324264098</li>\n",
       "</ol>\n",
       "</dd>\n",
       "\t<dt>$comparisons</dt>\n",
       "\t\t<dd><ol class=list-inline>\n",
       "\t<li>'M_Abx - M_Abx_FMT'</li>\n",
       "\t<li>'M_Abx - M_Control'</li>\n",
       "\t<li>'M_Abx_FMT - M_Control'</li>\n",
       "</ol>\n",
       "</dd>\n",
       "</dl>\n"
      ],
      "text/latex": [
       "\\begin{description}\n",
       "\\item[\\$chi2] 2.98748699414051\n",
       "\\item[\\$Z] \\begin{enumerate*}\n",
       "\\item -0.981494767538991\n",
       "\\item -1.72053810052189\n",
       "\\item -0.693154047118709\n",
       "\\end{enumerate*}\n",
       "\n",
       "\\item[\\$P] \\begin{enumerate*}\n",
       "\\item 0.163174406730601\n",
       "\\item 0.0426673370880744\n",
       "\\item 0.244106441421366\n",
       "\\end{enumerate*}\n",
       "\n",
       "\\item[\\$P.adjusted] \\begin{enumerate*}\n",
       "\\item 0.489523220191803\n",
       "\\item 0.128002011264223\n",
       "\\item 0.732319324264098\n",
       "\\end{enumerate*}\n",
       "\n",
       "\\item[\\$comparisons] \\begin{enumerate*}\n",
       "\\item 'M\\_Abx - M\\_Abx\\_FMT'\n",
       "\\item 'M\\_Abx - M\\_Control'\n",
       "\\item 'M\\_Abx\\_FMT - M\\_Control'\n",
       "\\end{enumerate*}\n",
       "\n",
       "\\end{description}\n"
      ],
      "text/markdown": [
       "$chi2\n",
       ":   2.98748699414051\n",
       "$Z\n",
       ":   1. -0.981494767538991\n",
       "2. -1.72053810052189\n",
       "3. -0.693154047118709\n",
       "\n",
       "\n",
       "\n",
       "$P\n",
       ":   1. 0.163174406730601\n",
       "2. 0.0426673370880744\n",
       "3. 0.244106441421366\n",
       "\n",
       "\n",
       "\n",
       "$P.adjusted\n",
       ":   1. 0.489523220191803\n",
       "2. 0.128002011264223\n",
       "3. 0.732319324264098\n",
       "\n",
       "\n",
       "\n",
       "$comparisons\n",
       ":   1. 'M_Abx - M_Abx_FMT'\n",
       "2. 'M_Abx - M_Control'\n",
       "3. 'M_Abx_FMT - M_Control'\n",
       "\n",
       "\n",
       "\n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "$chi2\n",
       "[1] 2.987487\n",
       "\n",
       "$Z\n",
       "[1] -0.9814948 -1.7205381 -0.6931540\n",
       "\n",
       "$P\n",
       "[1] 0.16317441 0.04266734 0.24410644\n",
       "\n",
       "$P.adjusted\n",
       "[1] 0.4895232 0.1280020 0.7323193\n",
       "\n",
       "$comparisons\n",
       "[1] \"M_Abx - M_Abx_FMT\"     \"M_Abx - M_Control\"     \"M_Abx_FMT - M_Control\"\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "#Observed richness by treatment group - MALES - CECAL\n",
    "\n",
    "obs = estimate_richness(na.omit(alz_nt_cec_M), split = TRUE, measures = \"Observed\")\n",
    "obs$SampleID = rownames(obs)\n",
    "\n",
    "dfobs = merge(obs, sample_tab, on=\"SampleID\")\n",
    "\n",
    "ordobs = obs %>% group_by(\"Treatment\") %>%\n",
    "    summarize(med.obs = median(Observed)) %>%\n",
    "    arrange(desc(med.obs)) %>%\n",
    "    magrittr::extract2(\"Treatment\")\n",
    "\n",
    "stats.obs = dfobs %>%\n",
    "    group_by(Treatment) %>%\n",
    "    summarize(median = median(Observed), \n",
    "              mean = mean(Observed),\n",
    "              sd = sd(Observed),\n",
    "              sem = sd(Observed/sqrt(length(Observed)))) #%>%\n",
    "    #arrange(-desc(mean))\n",
    "stats.obs\n",
    "\n",
    "#CHECK FOR SIGNIFICANCE\n",
    "stats.krusk = kruskal.test(Observed~Treatment, data=dfobs)\n",
    "stats.krusk\n",
    "\n",
    "#POST-HOC ANALYSIS\n",
    "dunn = dunn.test(dfobs$Observed, g=dfobs$Treatment, method=\"bonferroni\")\n",
    "dunn"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 92,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th scope=col>Treatment</th><th scope=col>median</th><th scope=col>mean</th><th scope=col>sd</th><th scope=col>sem</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><td>M_Abx     </td><td>3.736119  </td><td>3.742036  </td><td>0.2182130 </td><td>0.06900502</td></tr>\n",
       "\t<tr><td>M_Abx_FMT </td><td>3.726856  </td><td>3.567437  </td><td>0.2657537 </td><td>0.08858457</td></tr>\n",
       "\t<tr><td>M_Control </td><td>3.668383  </td><td>3.800824  </td><td>0.4051074 </td><td>0.12810620</td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|lllll}\n",
       " Treatment & median & mean & sd & sem\\\\\n",
       "\\hline\n",
       "\t M\\_Abx     & 3.736119    & 3.742036    & 0.2182130   & 0.06900502 \\\\\n",
       "\t M\\_Abx\\_FMT & 3.726856      & 3.567437      & 0.2657537     & 0.08858457   \\\\\n",
       "\t M\\_Control & 3.668383    & 3.800824    & 0.4051074   & 0.12810620 \\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "| Treatment | median | mean | sd | sem |\n",
       "|---|---|---|---|---|\n",
       "| M_Abx      | 3.736119   | 3.742036   | 0.2182130  | 0.06900502 |\n",
       "| M_Abx_FMT  | 3.726856   | 3.567437   | 0.2657537  | 0.08858457 |\n",
       "| M_Control  | 3.668383   | 3.800824   | 0.4051074  | 0.12810620 |\n",
       "\n"
      ],
      "text/plain": [
       "  Treatment median   mean     sd        sem       \n",
       "1 M_Abx     3.736119 3.742036 0.2182130 0.06900502\n",
       "2 M_Abx_FMT 3.726856 3.567437 0.2657537 0.08858457\n",
       "3 M_Control 3.668383 3.800824 0.4051074 0.12810620"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "\n",
       "\tKruskal-Wallis rank sum test\n",
       "\n",
       "data:  Shannon by Treatment\n",
       "Kruskal-Wallis chi-squared = 1.9027, df = 2, p-value = 0.3862\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  Kruskal-Wallis rank sum test\n",
      "\n",
      "data: x and group\n",
      "Kruskal-Wallis chi-squared = 1.9027, df = 2, p-value = 0.39\n",
      "\n",
      "\n",
      "                           Comparison of x by group                            \n",
      "                                 (Bonferroni)                                  \n",
      "Col Mean-|\n",
      "Row Mean |      M_Abx   M_Abx_FM\n",
      "---------+----------------------\n",
      "M_Abx_FM |   1.283723\n",
      "         |     0.2989\n",
      "         |\n",
      "M_Contro |   0.183829  -1.104797\n",
      "         |     1.0000     0.4039\n",
      "\n",
      "alpha = 0.05\n",
      "Reject Ho if p <= alpha/2\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<dl>\n",
       "\t<dt>$chi2</dt>\n",
       "\t\t<dd>1.90268199233716</dd>\n",
       "\t<dt>$Z</dt>\n",
       "\t\t<dd><ol class=list-inline>\n",
       "\t<li>1.28372315949815</li>\n",
       "\t<li>0.183829006003611</li>\n",
       "\t<li>-1.10479714390438</li>\n",
       "</ol>\n",
       "</dd>\n",
       "\t<dt>$P</dt>\n",
       "\t\t<dd><ol class=list-inline>\n",
       "\t<li>0.0996194190087712</li>\n",
       "\t<li>0.427073800196915</li>\n",
       "\t<li>0.134623749190596</li>\n",
       "</ol>\n",
       "</dd>\n",
       "\t<dt>$P.adjusted</dt>\n",
       "\t\t<dd><ol class=list-inline>\n",
       "\t<li>0.298858257026314</li>\n",
       "\t<li>1</li>\n",
       "\t<li>0.403871247571789</li>\n",
       "</ol>\n",
       "</dd>\n",
       "\t<dt>$comparisons</dt>\n",
       "\t\t<dd><ol class=list-inline>\n",
       "\t<li>'M_Abx - M_Abx_FMT'</li>\n",
       "\t<li>'M_Abx - M_Control'</li>\n",
       "\t<li>'M_Abx_FMT - M_Control'</li>\n",
       "</ol>\n",
       "</dd>\n",
       "</dl>\n"
      ],
      "text/latex": [
       "\\begin{description}\n",
       "\\item[\\$chi2] 1.90268199233716\n",
       "\\item[\\$Z] \\begin{enumerate*}\n",
       "\\item 1.28372315949815\n",
       "\\item 0.183829006003611\n",
       "\\item -1.10479714390438\n",
       "\\end{enumerate*}\n",
       "\n",
       "\\item[\\$P] \\begin{enumerate*}\n",
       "\\item 0.0996194190087712\n",
       "\\item 0.427073800196915\n",
       "\\item 0.134623749190596\n",
       "\\end{enumerate*}\n",
       "\n",
       "\\item[\\$P.adjusted] \\begin{enumerate*}\n",
       "\\item 0.298858257026314\n",
       "\\item 1\n",
       "\\item 0.403871247571789\n",
       "\\end{enumerate*}\n",
       "\n",
       "\\item[\\$comparisons] \\begin{enumerate*}\n",
       "\\item 'M\\_Abx - M\\_Abx\\_FMT'\n",
       "\\item 'M\\_Abx - M\\_Control'\n",
       "\\item 'M\\_Abx\\_FMT - M\\_Control'\n",
       "\\end{enumerate*}\n",
       "\n",
       "\\end{description}\n"
      ],
      "text/markdown": [
       "$chi2\n",
       ":   1.90268199233716\n",
       "$Z\n",
       ":   1. 1.28372315949815\n",
       "2. 0.183829006003611\n",
       "3. -1.10479714390438\n",
       "\n",
       "\n",
       "\n",
       "$P\n",
       ":   1. 0.0996194190087712\n",
       "2. 0.427073800196915\n",
       "3. 0.134623749190596\n",
       "\n",
       "\n",
       "\n",
       "$P.adjusted\n",
       ":   1. 0.298858257026314\n",
       "2. 1\n",
       "3. 0.403871247571789\n",
       "\n",
       "\n",
       "\n",
       "$comparisons\n",
       ":   1. 'M_Abx - M_Abx_FMT'\n",
       "2. 'M_Abx - M_Control'\n",
       "3. 'M_Abx_FMT - M_Control'\n",
       "\n",
       "\n",
       "\n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "$chi2\n",
       "[1] 1.902682\n",
       "\n",
       "$Z\n",
       "[1]  1.283723  0.183829 -1.104797\n",
       "\n",
       "$P\n",
       "[1] 0.09961942 0.42707380 0.13462375\n",
       "\n",
       "$P.adjusted\n",
       "[1] 0.2988583 1.0000000 0.4038712\n",
       "\n",
       "$comparisons\n",
       "[1] \"M_Abx - M_Abx_FMT\"     \"M_Abx - M_Control\"     \"M_Abx_FMT - M_Control\"\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "#Shanon Index by treatment group - MALES - FECAL\n",
    "\n",
    "sha = estimate_richness(na.omit(alz_nt_fec_M), split = TRUE, measures = \"Shannon\")\n",
    "sha$SampleID = rownames(sha)\n",
    "\n",
    "dfsha = merge(sha, sample_tab, on=\"SampleID\")\n",
    "\n",
    "ordsha = sha %>% group_by(\"Treatment\") %>%\n",
    "    summarize(med.sha = median(Shannon)) %>%\n",
    "    arrange(desc(med.sha)) %>%\n",
    "    magrittr::extract2(\"Treatment\")\n",
    "\n",
    "stats.sha = dfsha %>%\n",
    "    group_by(Treatment) %>%\n",
    "    summarize(median = median(Shannon), \n",
    "              mean = mean(Shannon),\n",
    "              sd = sd(Shannon),\n",
    "              sem = sd(Shannon/sqrt(length(Shannon)))) #%>%\n",
    "    #arrange(-desc(mean))\n",
    "stats.sha\n",
    "\n",
    "#CHECK FOR SIGNIFICANCE\n",
    "stats.krusk = kruskal.test(Shannon~Treatment, data=dfsha)\n",
    "stats.krusk\n",
    "\n",
    "#POST-HOC ANALYSIS\n",
    "dunn = dunn.test(dfsha$Shannon, g=dfsha$Treatment, method=\"bonferroni\")\n",
    "dunn"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 93,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th scope=col>Treatment</th><th scope=col>median</th><th scope=col>mean</th><th scope=col>sd</th><th scope=col>sem</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><td>M_Abx     </td><td>3.190896  </td><td>3.259377  </td><td>0.2855643 </td><td>0.09030336</td></tr>\n",
       "\t<tr><td>M_Abx_FMT </td><td>3.146525  </td><td>3.229104  </td><td>0.6027147 </td><td>0.20090489</td></tr>\n",
       "\t<tr><td>M_Control </td><td>3.505665  </td><td>3.504035  </td><td>0.5817431 </td><td>0.18396332</td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|lllll}\n",
       " Treatment & median & mean & sd & sem\\\\\n",
       "\\hline\n",
       "\t M\\_Abx     & 3.190896    & 3.259377    & 0.2855643   & 0.09030336 \\\\\n",
       "\t M\\_Abx\\_FMT & 3.146525      & 3.229104      & 0.6027147     & 0.20090489   \\\\\n",
       "\t M\\_Control & 3.505665    & 3.504035    & 0.5817431   & 0.18396332 \\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "| Treatment | median | mean | sd | sem |\n",
       "|---|---|---|---|---|\n",
       "| M_Abx      | 3.190896   | 3.259377   | 0.2855643  | 0.09030336 |\n",
       "| M_Abx_FMT  | 3.146525   | 3.229104   | 0.6027147  | 0.20090489 |\n",
       "| M_Control  | 3.505665   | 3.504035   | 0.5817431  | 0.18396332 |\n",
       "\n"
      ],
      "text/plain": [
       "  Treatment median   mean     sd        sem       \n",
       "1 M_Abx     3.190896 3.259377 0.2855643 0.09030336\n",
       "2 M_Abx_FMT 3.146525 3.229104 0.6027147 0.20090489\n",
       "3 M_Control 3.505665 3.504035 0.5817431 0.18396332"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "\n",
       "\tKruskal-Wallis rank sum test\n",
       "\n",
       "data:  Shannon by Treatment\n",
       "Kruskal-Wallis chi-squared = 1.3661, df = 2, p-value = 0.5051\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  Kruskal-Wallis rank sum test\n",
      "\n",
      "data: x and group\n",
      "Kruskal-Wallis chi-squared = 1.3661, df = 2, p-value = 0.51\n",
      "\n",
      "\n",
      "                           Comparison of x by group                            \n",
      "                                 (Bonferroni)                                  \n",
      "Col Mean-|\n",
      "Row Mean |      M_Abx   M_Abx_FM\n",
      "---------+----------------------\n",
      "M_Abx_FM |   0.224367\n",
      "         |     1.0000\n",
      "         |\n",
      "M_Contro |  -0.892883  -1.093436\n",
      "         |     0.5579     0.4113\n",
      "\n",
      "alpha = 0.05\n",
      "Reject Ho if p <= alpha/2\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<dl>\n",
       "\t<dt>$chi2</dt>\n",
       "\t\t<dd>1.36613026819923</dd>\n",
       "\t<dt>$Z</dt>\n",
       "\t\t<dd><ol class=list-inline>\n",
       "\t<li>0.224367543363614</li>\n",
       "\t<li>-0.892883743446113</li>\n",
       "\t<li>-1.09343676196192</li>\n",
       "</ol>\n",
       "</dd>\n",
       "\t<dt>$P</dt>\n",
       "\t\t<dd><ol class=list-inline>\n",
       "\t<li>0.411235661775699</li>\n",
       "\t<li>0.185959716840888</li>\n",
       "\t<li>0.137101040085233</li>\n",
       "</ol>\n",
       "</dd>\n",
       "\t<dt>$P.adjusted</dt>\n",
       "\t\t<dd><ol class=list-inline>\n",
       "\t<li>1</li>\n",
       "\t<li>0.557879150522665</li>\n",
       "\t<li>0.411303120255699</li>\n",
       "</ol>\n",
       "</dd>\n",
       "\t<dt>$comparisons</dt>\n",
       "\t\t<dd><ol class=list-inline>\n",
       "\t<li>'M_Abx - M_Abx_FMT'</li>\n",
       "\t<li>'M_Abx - M_Control'</li>\n",
       "\t<li>'M_Abx_FMT - M_Control'</li>\n",
       "</ol>\n",
       "</dd>\n",
       "</dl>\n"
      ],
      "text/latex": [
       "\\begin{description}\n",
       "\\item[\\$chi2] 1.36613026819923\n",
       "\\item[\\$Z] \\begin{enumerate*}\n",
       "\\item 0.224367543363614\n",
       "\\item -0.892883743446113\n",
       "\\item -1.09343676196192\n",
       "\\end{enumerate*}\n",
       "\n",
       "\\item[\\$P] \\begin{enumerate*}\n",
       "\\item 0.411235661775699\n",
       "\\item 0.185959716840888\n",
       "\\item 0.137101040085233\n",
       "\\end{enumerate*}\n",
       "\n",
       "\\item[\\$P.adjusted] \\begin{enumerate*}\n",
       "\\item 1\n",
       "\\item 0.557879150522665\n",
       "\\item 0.411303120255699\n",
       "\\end{enumerate*}\n",
       "\n",
       "\\item[\\$comparisons] \\begin{enumerate*}\n",
       "\\item 'M\\_Abx - M\\_Abx\\_FMT'\n",
       "\\item 'M\\_Abx - M\\_Control'\n",
       "\\item 'M\\_Abx\\_FMT - M\\_Control'\n",
       "\\end{enumerate*}\n",
       "\n",
       "\\end{description}\n"
      ],
      "text/markdown": [
       "$chi2\n",
       ":   1.36613026819923\n",
       "$Z\n",
       ":   1. 0.224367543363614\n",
       "2. -0.892883743446113\n",
       "3. -1.09343676196192\n",
       "\n",
       "\n",
       "\n",
       "$P\n",
       ":   1. 0.411235661775699\n",
       "2. 0.185959716840888\n",
       "3. 0.137101040085233\n",
       "\n",
       "\n",
       "\n",
       "$P.adjusted\n",
       ":   1. 1\n",
       "2. 0.557879150522665\n",
       "3. 0.411303120255699\n",
       "\n",
       "\n",
       "\n",
       "$comparisons\n",
       ":   1. 'M_Abx - M_Abx_FMT'\n",
       "2. 'M_Abx - M_Control'\n",
       "3. 'M_Abx_FMT - M_Control'\n",
       "\n",
       "\n",
       "\n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "$chi2\n",
       "[1] 1.36613\n",
       "\n",
       "$Z\n",
       "[1]  0.2243675 -0.8928837 -1.0934368\n",
       "\n",
       "$P\n",
       "[1] 0.4112357 0.1859597 0.1371010\n",
       "\n",
       "$P.adjusted\n",
       "[1] 1.0000000 0.5578792 0.4113031\n",
       "\n",
       "$comparisons\n",
       "[1] \"M_Abx - M_Abx_FMT\"     \"M_Abx - M_Control\"     \"M_Abx_FMT - M_Control\"\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "#Shanon Index by treatment group - MALES - CECAL\n",
    "\n",
    "sha = estimate_richness(na.omit(alz_nt_cec_M), split = TRUE, measures = \"Shannon\")\n",
    "sha$SampleID = rownames(sha)\n",
    "\n",
    "dfsha = merge(sha, sample_tab, on=\"SampleID\")\n",
    "\n",
    "ordsha = sha %>% group_by(\"Treatment\") %>%\n",
    "    summarize(med.sha = median(Shannon)) %>%\n",
    "    arrange(desc(med.sha)) %>%\n",
    "    magrittr::extract2(\"Treatment\")\n",
    "\n",
    "stats.sha = dfsha %>%\n",
    "    group_by(Treatment) %>%\n",
    "    summarize(median = median(Shannon), \n",
    "              mean = mean(Shannon),\n",
    "              sd = sd(Shannon),\n",
    "              sem = sd(Shannon/sqrt(length(Shannon)))) #%>%\n",
    "    #arrange(-desc(mean))\n",
    "stats.sha\n",
    "\n",
    "#CHECK FOR SIGNIFICANCE\n",
    "stats.krusk = kruskal.test(Shannon~Treatment, data=dfsha)\n",
    "stats.krusk\n",
    "\n",
    "#POST-HOC ANALYSIS\n",
    "dunn = dunn.test(dfsha$Shannon, g=dfsha$Treatment, method=\"bonferroni\")\n",
    "dunn"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 90,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th scope=col>Treatment</th><th scope=col>median</th><th scope=col>mean</th><th scope=col>sd</th><th scope=col>sem</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><td>F_Abx    </td><td>141      </td><td>143.6667 </td><td>40.42586 </td><td>13.47529 </td></tr>\n",
       "\t<tr><td>F_Control</td><td>203      </td><td>199.9000 </td><td>41.60248 </td><td>13.15586 </td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|lllll}\n",
       " Treatment & median & mean & sd & sem\\\\\n",
       "\\hline\n",
       "\t F\\_Abx     & 141         & 143.6667    & 40.42586    & 13.47529   \\\\\n",
       "\t F\\_Control & 203         & 199.9000    & 41.60248    & 13.15586   \\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "| Treatment | median | mean | sd | sem |\n",
       "|---|---|---|---|---|\n",
       "| F_Abx     | 141       | 143.6667  | 40.42586  | 13.47529  |\n",
       "| F_Control | 203       | 199.9000  | 41.60248  | 13.15586  |\n",
       "\n"
      ],
      "text/plain": [
       "  Treatment median mean     sd       sem     \n",
       "1 F_Abx     141    143.6667 40.42586 13.47529\n",
       "2 F_Control 203    199.9000 41.60248 13.15586"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "\n",
       "\tKruskal-Wallis rank sum test\n",
       "\n",
       "data:  Observed by Treatment\n",
       "Kruskal-Wallis chi-squared = 6.4123, df = 1, p-value = 0.01133\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  Kruskal-Wallis rank sum test\n",
      "\n",
      "data: x and group\n",
      "Kruskal-Wallis chi-squared = 6.4123, df = 1, p-value = 0.01\n",
      "\n",
      "\n",
      "                           Comparison of x by group                            \n",
      "                                 (Bonferroni)                                  \n",
      "Col Mean-|\n",
      "Row Mean |      F_Abx\n",
      "---------+-----------\n",
      "F_Contro |  -2.532250\n",
      "         |    0.0057*\n",
      "\n",
      "alpha = 0.05\n",
      "Reject Ho if p <= alpha/2\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<dl>\n",
       "\t<dt>$chi2</dt>\n",
       "\t\t<dd>6.41229148375768</dd>\n",
       "\t<dt>$Z</dt>\n",
       "\t\t<dd>-2.53225028063138</dd>\n",
       "\t<dt>$P</dt>\n",
       "\t\t<dd>0.00566665321181593</dd>\n",
       "\t<dt>$P.adjusted</dt>\n",
       "\t\t<dd>0.00566665321181593</dd>\n",
       "\t<dt>$comparisons</dt>\n",
       "\t\t<dd>'F_Abx - F_Control'</dd>\n",
       "</dl>\n"
      ],
      "text/latex": [
       "\\begin{description}\n",
       "\\item[\\$chi2] 6.41229148375768\n",
       "\\item[\\$Z] -2.53225028063138\n",
       "\\item[\\$P] 0.00566665321181593\n",
       "\\item[\\$P.adjusted] 0.00566665321181593\n",
       "\\item[\\$comparisons] 'F\\_Abx - F\\_Control'\n",
       "\\end{description}\n"
      ],
      "text/markdown": [
       "$chi2\n",
       ":   6.41229148375768\n",
       "$Z\n",
       ":   -2.53225028063138\n",
       "$P\n",
       ":   0.00566665321181593\n",
       "$P.adjusted\n",
       ":   0.00566665321181593\n",
       "$comparisons\n",
       ":   'F_Abx - F_Control'\n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "$chi2\n",
       "[1] 6.412291\n",
       "\n",
       "$Z\n",
       "[1] -2.53225\n",
       "\n",
       "$P\n",
       "[1] 0.005666653\n",
       "\n",
       "$P.adjusted\n",
       "[1] 0.005666653\n",
       "\n",
       "$comparisons\n",
       "[1] \"F_Abx - F_Control\"\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "#Observed richness by treatment group - FEMALES - FECAL\n",
    "\n",
    "obs = estimate_richness(na.omit(alz_nt_fec_F), split = TRUE, measures = \"Observed\")\n",
    "obs$SampleID = rownames(obs)\n",
    "\n",
    "dfobs = merge(obs, sample_tab, on=\"SampleID\")\n",
    "\n",
    "ordobs = obs %>% group_by(\"Treatment\") %>%\n",
    "    summarize(med.obs = median(Observed)) %>%\n",
    "    arrange(desc(med.obs)) %>%\n",
    "    magrittr::extract2(\"Treatment\")\n",
    "\n",
    "stats.obs = dfobs %>%\n",
    "    group_by(Treatment) %>%\n",
    "    summarize(median = median(Observed), \n",
    "              mean = mean(Observed),\n",
    "              sd = sd(Observed),\n",
    "              sem = sd(Observed/sqrt(length(Observed)))) #%>%\n",
    "    #arrange(-desc(mean))\n",
    "stats.obs\n",
    "\n",
    "#CHECK FOR SIGNIFICANCE\n",
    "stats.krusk = kruskal.test(Observed~Treatment, data=dfobs)\n",
    "stats.krusk\n",
    "\n",
    "#POST-HOC ANALYSIS\n",
    "dunn = dunn.test(dfobs$Observed, g=dfobs$Treatment, method=\"bonferroni\")\n",
    "dunn"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 91,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th scope=col>Treatment</th><th scope=col>median</th><th scope=col>mean</th><th scope=col>sd</th><th scope=col>sem</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><td>F_Abx    </td><td>108.0    </td><td>118.875  </td><td>43.66901 </td><td>15.43933 </td></tr>\n",
       "\t<tr><td>F_Control</td><td>229.5    </td><td>230.500  </td><td>37.21783 </td><td>11.76931 </td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|lllll}\n",
       " Treatment & median & mean & sd & sem\\\\\n",
       "\\hline\n",
       "\t F\\_Abx     & 108.0       & 118.875     & 43.66901    & 15.43933   \\\\\n",
       "\t F\\_Control & 229.5       & 230.500     & 37.21783    & 11.76931   \\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "| Treatment | median | mean | sd | sem |\n",
       "|---|---|---|---|---|\n",
       "| F_Abx     | 108.0     | 118.875   | 43.66901  | 15.43933  |\n",
       "| F_Control | 229.5     | 230.500   | 37.21783  | 11.76931  |\n",
       "\n"
      ],
      "text/plain": [
       "  Treatment median mean    sd       sem     \n",
       "1 F_Abx     108.0  118.875 43.66901 15.43933\n",
       "2 F_Control 229.5  230.500 37.21783 11.76931"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "\n",
       "\tKruskal-Wallis rank sum test\n",
       "\n",
       "data:  Observed by Treatment\n",
       "Kruskal-Wallis chi-squared = 10.853, df = 1, p-value = 0.0009865\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  Kruskal-Wallis rank sum test\n",
      "\n",
      "data: x and group\n",
      "Kruskal-Wallis chi-squared = 10.8527, df = 1, p-value = 0\n",
      "\n",
      "\n",
      "                           Comparison of x by group                            \n",
      "                                 (Bonferroni)                                  \n",
      "Col Mean-|\n",
      "Row Mean |      F_Abx\n",
      "---------+-----------\n",
      "F_Contro |  -3.294342\n",
      "         |    0.0005*\n",
      "\n",
      "alpha = 0.05\n",
      "Reject Ho if p <= alpha/2\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<dl>\n",
       "\t<dt>$chi2</dt>\n",
       "\t\t<dd>10.8526943005181</dd>\n",
       "\t<dt>$Z</dt>\n",
       "\t\t<dd>-3.29434277216536</dd>\n",
       "\t<dt>$P</dt>\n",
       "\t\t<dd>0.000493260586987488</dd>\n",
       "\t<dt>$P.adjusted</dt>\n",
       "\t\t<dd>0.000493260586987488</dd>\n",
       "\t<dt>$comparisons</dt>\n",
       "\t\t<dd>'F_Abx - F_Control'</dd>\n",
       "</dl>\n"
      ],
      "text/latex": [
       "\\begin{description}\n",
       "\\item[\\$chi2] 10.8526943005181\n",
       "\\item[\\$Z] -3.29434277216536\n",
       "\\item[\\$P] 0.000493260586987488\n",
       "\\item[\\$P.adjusted] 0.000493260586987488\n",
       "\\item[\\$comparisons] 'F\\_Abx - F\\_Control'\n",
       "\\end{description}\n"
      ],
      "text/markdown": [
       "$chi2\n",
       ":   10.8526943005181\n",
       "$Z\n",
       ":   -3.29434277216536\n",
       "$P\n",
       ":   0.000493260586987488\n",
       "$P.adjusted\n",
       ":   0.000493260586987488\n",
       "$comparisons\n",
       ":   'F_Abx - F_Control'\n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "$chi2\n",
       "[1] 10.85269\n",
       "\n",
       "$Z\n",
       "[1] -3.294343\n",
       "\n",
       "$P\n",
       "[1] 0.0004932606\n",
       "\n",
       "$P.adjusted\n",
       "[1] 0.0004932606\n",
       "\n",
       "$comparisons\n",
       "[1] \"F_Abx - F_Control\"\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "#Observed richness by treatment group - FEMALES - CECAL\n",
    "\n",
    "obs = estimate_richness(na.omit(alz_nt_cec_F), split = TRUE, measures = \"Observed\")\n",
    "obs$SampleID = rownames(obs)\n",
    "\n",
    "dfobs = merge(obs, sample_tab, on=\"SampleID\")\n",
    "\n",
    "ordobs = obs %>% group_by(\"Treatment\") %>%\n",
    "    summarize(med.obs = median(Observed)) %>%\n",
    "    arrange(desc(med.obs)) %>%\n",
    "    magrittr::extract2(\"Treatment\")\n",
    "\n",
    "stats.obs = dfobs %>%\n",
    "    group_by(Treatment) %>%\n",
    "    summarize(median = median(Observed), \n",
    "              mean = mean(Observed),\n",
    "              sd = sd(Observed),\n",
    "              sem = sd(Observed/sqrt(length(Observed)))) #%>%\n",
    "    #arrange(-desc(mean))\n",
    "stats.obs\n",
    "\n",
    "#CHECK FOR SIGNIFICANCE\n",
    "stats.krusk = kruskal.test(Observed~Treatment, data=dfobs)\n",
    "stats.krusk\n",
    "\n",
    "#POST-HOC ANALYSIS\n",
    "dunn = dunn.test(dfobs$Observed, g=dfobs$Treatment, method=\"bonferroni\")\n",
    "dunn"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 94,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th scope=col>Treatment</th><th scope=col>median</th><th scope=col>mean</th><th scope=col>sd</th><th scope=col>sem</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><td>F_Abx    </td><td>3.524927 </td><td>3.671012 </td><td>0.3272313</td><td>0.1090771</td></tr>\n",
       "\t<tr><td>F_Control</td><td>3.762176 </td><td>3.864489 </td><td>0.3151918</td><td>0.0996724</td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|lllll}\n",
       " Treatment & median & mean & sd & sem\\\\\n",
       "\\hline\n",
       "\t F\\_Abx     & 3.524927    & 3.671012    & 0.3272313   & 0.1090771  \\\\\n",
       "\t F\\_Control & 3.762176    & 3.864489    & 0.3151918   & 0.0996724  \\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "| Treatment | median | mean | sd | sem |\n",
       "|---|---|---|---|---|\n",
       "| F_Abx     | 3.524927  | 3.671012  | 0.3272313 | 0.1090771 |\n",
       "| F_Control | 3.762176  | 3.864489  | 0.3151918 | 0.0996724 |\n",
       "\n"
      ],
      "text/plain": [
       "  Treatment median   mean     sd        sem      \n",
       "1 F_Abx     3.524927 3.671012 0.3272313 0.1090771\n",
       "2 F_Control 3.762176 3.864489 0.3151918 0.0996724"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "\n",
       "\tKruskal-Wallis rank sum test\n",
       "\n",
       "data:  Shannon by Treatment\n",
       "Kruskal-Wallis chi-squared = 1.9267, df = 1, p-value = 0.1651\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  Kruskal-Wallis rank sum test\n",
      "\n",
      "data: x and group\n",
      "Kruskal-Wallis chi-squared = 1.9267, df = 1, p-value = 0.17\n",
      "\n",
      "\n",
      "                           Comparison of x by group                            \n",
      "                                 (Bonferroni)                                  \n",
      "Col Mean-|\n",
      "Row Mean |      F_Abx\n",
      "---------+-----------\n",
      "F_Contro |  -1.388044\n",
      "         |     0.0826\n",
      "\n",
      "alpha = 0.05\n",
      "Reject Ho if p <= alpha/2\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<dl>\n",
       "\t<dt>$chi2</dt>\n",
       "\t\t<dd>1.92666666666667</dd>\n",
       "\t<dt>$Z</dt>\n",
       "\t\t<dd>-1.38804418757713</dd>\n",
       "\t<dt>$P</dt>\n",
       "\t\t<dd>0.0825617950641971</dd>\n",
       "\t<dt>$P.adjusted</dt>\n",
       "\t\t<dd>0.0825617950641971</dd>\n",
       "\t<dt>$comparisons</dt>\n",
       "\t\t<dd>'F_Abx - F_Control'</dd>\n",
       "</dl>\n"
      ],
      "text/latex": [
       "\\begin{description}\n",
       "\\item[\\$chi2] 1.92666666666667\n",
       "\\item[\\$Z] -1.38804418757713\n",
       "\\item[\\$P] 0.0825617950641971\n",
       "\\item[\\$P.adjusted] 0.0825617950641971\n",
       "\\item[\\$comparisons] 'F\\_Abx - F\\_Control'\n",
       "\\end{description}\n"
      ],
      "text/markdown": [
       "$chi2\n",
       ":   1.92666666666667\n",
       "$Z\n",
       ":   -1.38804418757713\n",
       "$P\n",
       ":   0.0825617950641971\n",
       "$P.adjusted\n",
       ":   0.0825617950641971\n",
       "$comparisons\n",
       ":   'F_Abx - F_Control'\n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "$chi2\n",
       "[1] 1.926667\n",
       "\n",
       "$Z\n",
       "[1] -1.388044\n",
       "\n",
       "$P\n",
       "[1] 0.0825618\n",
       "\n",
       "$P.adjusted\n",
       "[1] 0.0825618\n",
       "\n",
       "$comparisons\n",
       "[1] \"F_Abx - F_Control\"\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "#Shanon Index by treatment group - FEMALES - FECAL\n",
    "\n",
    "sha = estimate_richness(na.omit(alz_nt_fec_F), split = TRUE, measures = \"Shannon\")\n",
    "sha$SampleID = rownames(sha)\n",
    "\n",
    "dfsha = merge(sha, sample_tab, on=\"SampleID\")\n",
    "\n",
    "ordsha = sha %>% group_by(\"Treatment\") %>%\n",
    "    summarize(med.sha = median(Shannon)) %>%\n",
    "    arrange(desc(med.sha)) %>%\n",
    "    magrittr::extract2(\"Treatment\")\n",
    "\n",
    "stats.sha = dfsha %>%\n",
    "    group_by(Treatment) %>%\n",
    "    summarize(median = median(Shannon), \n",
    "              mean = mean(Shannon),\n",
    "              sd = sd(Shannon),\n",
    "              sem = sd(Shannon/sqrt(length(Shannon)))) #%>%\n",
    "    #arrange(-desc(mean))\n",
    "stats.sha\n",
    "\n",
    "#CHECK FOR SIGNIFICANCE\n",
    "stats.krusk = kruskal.test(Shannon~Treatment, data=dfsha)\n",
    "stats.krusk\n",
    "\n",
    "#POST-HOC ANALYSIS\n",
    "dunn = dunn.test(dfsha$Shannon, g=dfsha$Treatment, method=\"bonferroni\")\n",
    "dunn"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 95,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th scope=col>Treatment</th><th scope=col>median</th><th scope=col>mean</th><th scope=col>sd</th><th scope=col>sem</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><td>F_Abx    </td><td>3.552092 </td><td>3.352493 </td><td>0.4709796</td><td>0.1665164</td></tr>\n",
       "\t<tr><td>F_Control</td><td>3.989772 </td><td>3.947161 </td><td>0.3168112</td><td>0.1001845</td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|lllll}\n",
       " Treatment & median & mean & sd & sem\\\\\n",
       "\\hline\n",
       "\t F\\_Abx     & 3.552092    & 3.352493    & 0.4709796   & 0.1665164  \\\\\n",
       "\t F\\_Control & 3.989772    & 3.947161    & 0.3168112   & 0.1001845  \\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "| Treatment | median | mean | sd | sem |\n",
       "|---|---|---|---|---|\n",
       "| F_Abx     | 3.552092  | 3.352493  | 0.4709796 | 0.1665164 |\n",
       "| F_Control | 3.989772  | 3.947161  | 0.3168112 | 0.1001845 |\n",
       "\n"
      ],
      "text/plain": [
       "  Treatment median   mean     sd        sem      \n",
       "1 F_Abx     3.552092 3.352493 0.4709796 0.1665164\n",
       "2 F_Control 3.989772 3.947161 0.3168112 0.1001845"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "\n",
       "\tKruskal-Wallis rank sum test\n",
       "\n",
       "data:  Shannon by Treatment\n",
       "Kruskal-Wallis chi-squared = 8.5974, df = 1, p-value = 0.003366\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  Kruskal-Wallis rank sum test\n",
      "\n",
      "data: x and group\n",
      "Kruskal-Wallis chi-squared = 8.5974, df = 1, p-value = 0\n",
      "\n",
      "\n",
      "                           Comparison of x by group                            \n",
      "                                 (Bonferroni)                                  \n",
      "Col Mean-|\n",
      "Row Mean |      F_Abx\n",
      "---------+-----------\n",
      "F_Contro |  -2.932126\n",
      "         |    0.0017*\n",
      "\n",
      "alpha = 0.05\n",
      "Reject Ho if p <= alpha/2\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<dl>\n",
       "\t<dt>$chi2</dt>\n",
       "\t\t<dd>8.59736842105264</dd>\n",
       "\t<dt>$Z</dt>\n",
       "\t\t<dd>-2.93212694490751</dd>\n",
       "\t<dt>$P</dt>\n",
       "\t\t<dd>0.0016832455399728</dd>\n",
       "\t<dt>$P.adjusted</dt>\n",
       "\t\t<dd>0.0016832455399728</dd>\n",
       "\t<dt>$comparisons</dt>\n",
       "\t\t<dd>'F_Abx - F_Control'</dd>\n",
       "</dl>\n"
      ],
      "text/latex": [
       "\\begin{description}\n",
       "\\item[\\$chi2] 8.59736842105264\n",
       "\\item[\\$Z] -2.93212694490751\n",
       "\\item[\\$P] 0.0016832455399728\n",
       "\\item[\\$P.adjusted] 0.0016832455399728\n",
       "\\item[\\$comparisons] 'F\\_Abx - F\\_Control'\n",
       "\\end{description}\n"
      ],
      "text/markdown": [
       "$chi2\n",
       ":   8.59736842105264\n",
       "$Z\n",
       ":   -2.93212694490751\n",
       "$P\n",
       ":   0.0016832455399728\n",
       "$P.adjusted\n",
       ":   0.0016832455399728\n",
       "$comparisons\n",
       ":   'F_Abx - F_Control'\n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "$chi2\n",
       "[1] 8.597368\n",
       "\n",
       "$Z\n",
       "[1] -2.932127\n",
       "\n",
       "$P\n",
       "[1] 0.001683246\n",
       "\n",
       "$P.adjusted\n",
       "[1] 0.001683246\n",
       "\n",
       "$comparisons\n",
       "[1] \"F_Abx - F_Control\"\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "#Shanon Index by treatment group - FEMALES - CECAL\n",
    "\n",
    "sha = estimate_richness(na.omit(alz_nt_cec_F), split = TRUE, measures = \"Shannon\")\n",
    "sha$SampleID = rownames(sha)\n",
    "\n",
    "dfsha = merge(sha, sample_tab, on=\"SampleID\")\n",
    "\n",
    "ordsha = sha %>% group_by(\"Treatment\") %>%\n",
    "    summarize(med.sha = median(Shannon)) %>%\n",
    "    arrange(desc(med.sha)) %>%\n",
    "    magrittr::extract2(\"Treatment\")\n",
    "\n",
    "stats.sha = dfsha %>%\n",
    "    group_by(Treatment) %>%\n",
    "    summarize(median = median(Shannon), \n",
    "              mean = mean(Shannon),\n",
    "              sd = sd(Shannon),\n",
    "              sem = sd(Shannon/sqrt(length(Shannon)))) #%>%\n",
    "    #arrange(-desc(mean))\n",
    "stats.sha\n",
    "\n",
    "#CHECK FOR SIGNIFICANCE\n",
    "stats.krusk = kruskal.test(Shannon~Treatment, data=dfsha)\n",
    "stats.krusk\n",
    "\n",
    "#POST-HOC ANALYSIS\n",
    "dunn = dunn.test(dfsha$Shannon, g=dfsha$Treatment, method=\"bonferroni\")\n",
    "dunn"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### UNWEIGHTED UNIFRAC"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [],
   "source": [
    "alz_t_M = subset_samples(alz_t_fec_, Sex==\"M\")\n",
    "alz_t_F = subset_samples(alz_t_1k, Sex==\"F\")\n",
    "\n",
    "alz_t_1k_M.df = data.frame(sample_data(alz_t_1k_M))\n",
    "alz_t_1k_F.df = data.frame(sample_data(alz_t_1k_F))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Warning message in UniFrac(physeq, ...):\n",
      "“Randomly assigning root as -- TACGTATGTCCCGAGCGTTATCCGGATTTATTGGGCGTAAAGCGAGCGCAGACGGTTGATTAAGTCTGATGTGAAAGCCCGGAGCTCAACCCCGGAATGGCATTGGAAACTGGTTAACTTGAGTGCAGTAGAGGTAAGTGGAACTCCATG -- in the phylogenetic tree in the data you provided.”"
     ]
    },
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th></th><th scope=col>Df</th><th scope=col>SumOfSqs</th><th scope=col>R2</th><th scope=col>F</th><th scope=col>Pr(&gt;F)</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><th scope=row>Treatment</th><td> 2       </td><td>0.4116223</td><td>0.0762443</td><td>1.114254 </td><td>0.193    </td></tr>\n",
       "\t<tr><th scope=row>Residual</th><td>27       </td><td>4.9871064</td><td>0.9237557</td><td>      NA </td><td>   NA    </td></tr>\n",
       "\t<tr><th scope=row>Total</th><td>29       </td><td>5.3987288</td><td>1.0000000</td><td>      NA </td><td>   NA    </td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|lllll}\n",
       "  & Df & SumOfSqs & R2 & F & Pr(>F)\\\\\n",
       "\\hline\n",
       "\tTreatment &  2        & 0.4116223 & 0.0762443 & 1.114254  & 0.193    \\\\\n",
       "\tResidual & 27        & 4.9871064 & 0.9237557 &       NA  &    NA    \\\\\n",
       "\tTotal & 29        & 5.3987288 & 1.0000000 &       NA  &    NA    \\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "| <!--/--> | Df | SumOfSqs | R2 | F | Pr(>F) |\n",
       "|---|---|---|---|---|---|\n",
       "| Treatment |  2        | 0.4116223 | 0.0762443 | 1.114254  | 0.193     |\n",
       "| Residual | 27        | 4.9871064 | 0.9237557 |       NA  |    NA     |\n",
       "| Total | 29        | 5.3987288 | 1.0000000 |       NA  |    NA     |\n",
       "\n"
      ],
      "text/plain": [
       "          Df SumOfSqs  R2        F        Pr(>F)\n",
       "Treatment  2 0.4116223 0.0762443 1.114254 0.193 \n",
       "Residual  27 4.9871064 0.9237557       NA    NA \n",
       "Total     29 5.3987288 1.0000000       NA    NA "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "alz.dist.unif = distance(alz_t_1k_M, method = \"unifrac\", type = \"samples\")\n",
    "alz.unif.treat.adtest = adonis2(alz.dist.unif ~ Treatment, data=alz_t_1k_M.df)\n",
    "alz.unif.treat.adtest"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "R_3.6.0",
   "language": "R",
   "name": "ir360"
  },
  "language_info": {
   "codemirror_mode": "r",
   "file_extension": ".r",
   "mimetype": "text/x-r-source",
   "name": "R",
   "pygments_lexer": "r",
   "version": "3.6.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}