class13_correlnet_R.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# CS446/546 - Class Session 19 - Correlation networks\n",
    "\n",
    "In this class session we are going to analyze gene expression data from a human bladder cancer cohort. We will load a data matrix of expression measurements of 4,473 genes in 414 different bladder cancer samples. These genes have been selected because they are differentially expressed between normal bladder and bladder cancer (thus more likely to have a function in bladder cancer specifically), but the columns in the data matrix are restricted to bladder cancer samples (not normal bladder) because we want to obtain a network representing variation across cancers. The measurements in the matrix have already been normalized to account for inter-sample heterogeneity and then log2 transformed. Our job is to compute Pearson correlation coefficients between all pairs of genes, obtain Fisher-transformed *z*-scores for all pairs of genes, test each pair of genes for significance of the z score, adjust for multiple hypothesis testing, filter to eliminate any pair for which *R* < 0.75 or *P*adj > 0.01, load the graph into an `igraph::Graph`, and plot the degree distribution on log-log scale. We will then answer two questions: (1) does the network look to be scale-free? and (2) what is it's best-fit scaling exponent?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's start by loading the packages that we will need for this notebook. Note the difference in language-design philosophy between R (which requires one package for this analysis) and python (where we have to load eight modules). Python keeps its core minimal, whereas R has a lot of statistical and plotting functions in the base language (or in packages that are loaded by default)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "suppressPackageStartupMessages(\n",
    "    library(igraph)\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Read the tab-deliminted text file of gene expression measurements (rows correspond to genes, columns correspond to bladder tumor samples). (use `read.table` with `row.names=1`). As always, sanity check that the file that you loaded has the expected dimensions (4,473 x 414)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<ol class=list-inline>\n",
       "\t<li>4473</li>\n",
       "\t<li>414</li>\n",
       "</ol>\n"
      ],
      "text/latex": [
       "\\begin{enumerate*}\n",
       "\\item 4473\n",
       "\\item 414\n",
       "\\end{enumerate*}\n"
      ],
      "text/markdown": [
       "1. 4473\n",
       "2. 414\n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "[1] 4473  414"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "gene_matrix_for_network <- read.table(\"shared/bladder_cancer_genes_tcga.txt\",\n",
    "                                     sep=\"\\t\",\n",
    "                                     header=TRUE,\n",
    "                                     row.names=1,\n",
    "                                     stringsAsFactors=FALSE)\n",
    "dim(gene_matrix_for_network)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Look up the online help for the `cor(x)` function, using the command `?cor`. When the `x` argument is a matrix, does `cor` compute the correlation coefficient of pairs of columns or pairs of rows?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "?cor"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since `cor(x)` computes the correlation coefficients of pairs of *columns* of its matrix argument but we want to compute correlation coefficients of pairs of *genes* (and recall that genes are *rows* in our matrix `gene_matrix_for_network`), use the `t` function to transpose the matrix. Then use the `cor` function to compute the 4,473 x 4,473 matrix of gene-gene Pearson correlation coefficients, and call this matrix `gene_matrix_for_network_cor`. It is in your interest to do this as a one-line command (by composing the function calls), so we store only one big matrix in the environment:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "gene_matrix_for_network_cor <- cor(t(gene_matrix_for_network))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Look up the online help for the `upper.tri` function. Does it return the matrix's upper triangle entries or the *index values* of the upper triangle entries of the matrix?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "?upper.tri"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Set the upper-triangle of the matrix `gene_matrix_for_network_cor` to zero, using (from left to right) three things: array indexing `[`, the function `upper.tri`, and the assignment operator `<-`. In the call to `upper.tri`, specify `diag=TRUE` to also zero out the diagonal."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "gene_matrix_for_network_cor[upper.tri(gene_matrix_for_network_cor, diag=TRUE)] <- 0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Use the function `which` in order to obtain a vector `inds_correl_above_thresh` containing the indices of the entries of the matrix `gene_matrix_for_network_cor` for which *R* >= 0.75. Use array indexing to obtain the *R* values for these matrix entries, as a numpy array `cor_coeff_values_above_thresh`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "inds_correl_above_thresh <- which(gene_matrix_for_network_cor >= 0.75)\n",
    "cor_coeff_values_above_thresh <- gene_matrix_for_network_cor[inds_correl_above_thresh]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Refer to Eq. (13.5) in the assigned readding for today's class (p9 of the PDF).  Obtain a vector of the correlation coefficients that exceeded 0.75, and Fisher-transform the correlation coefficient values to get a vector `z_scores` of *z* scores. Each of these *z* scores will correspond to an **edge** in the network, unless the absolute *z* score is too small such that we can't exclude the null hypothesis that the corresponding two genes' expression values are indepdenent (we will perform that check in the next step)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "z_scores <- 0.5*log((1+gene_matrix_for_network_cor[inds_correl_above_thresh])/\n",
    "                    (1-gene_matrix_for_network_cor[inds_correl_above_thresh]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Delete the correlation matrix object in order to save memory (we won't need it from here on out)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "rm(gene_matrix_for_network_cor)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we are going to perform a statistical test on *each* correlation coefficient that was >= 0.75.  Assume that under the null hypothesis that two genes are independent, then sqrt(M-3)z for the pair of genes is independent sample from the normal distribution with zero mean and unit variance, where M is the number of samples used to compute the Pearson correlation coefficient (i.e., M = 414). For each entry in `z_scores` compute a P value as the area under two tails of the normal distribution N(x),  where the two tails are x < -sqrt(M-3)z and x > sqrt(M-3)z. (You'll know you are doing it right if z=0 means you get a P value of 1).  You will want to use the functions `sqrt` and `pnorm`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "M <- ncol(gene_matrix_for_network)\n",
    "P_values <- 2*(pnorm(-z_scores*sqrt(M-3)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Look up the function `p.adjust` in the online help"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "?p.adjust"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Adjust the P values for multiple hypothesis testing, using the `p.adjust` function with `method=\"hochberg\"`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "P_values_adj <- p.adjust(P_values, method=\"hochberg\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Verify that we don't need to drop any entries due to the adjusted P value not being small enough (use `which` and `length`); this should produce zero since we have M=414 samples per gene."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "0"
      ],
      "text/latex": [
       "0"
      ],
      "text/markdown": [
       "0"
      ],
      "text/plain": [
       "[1] 0"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "length(which(P_values_adj > 0.01))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Make an undirected graph from the row/column indices of the (upper-triangle) gene pairs whose  correlations were above our threshold (you'll want to use the `arrayInd`, `cbind`, and `graph_from_data_frame` functions for this). Print a summary of the network, as a sanity check. Make sure to specify `directed=FALSE` when you construct the graph."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "IGRAPH b11fcc8 UN-- 1394 9916 -- \n",
      "+ attr: name (v/c)\n"
     ]
    }
   ],
   "source": [
    "N <- nrow(gene_matrix_for_network)\n",
    "final_network <- graph_from_data_frame(cbind(arrayInd(inds_correl_above_thresh, \n",
    "                                                      .dim=c(N,N))), directed=FALSE)\n",
    "summary(final_network)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Plot the degree distribution on log-log scale; does it appear to be scale-free?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAeAAAAHgCAMAAABKCk6nAAADAFBMVEUAAAABAQECAgIDAwME\nBAQFBQUGBgYHBwcICAgJCQkKCgoLCwsMDAwNDQ0ODg4PDw8QEBARERESEhITExMUFBQVFRUW\nFhYXFxcYGBgZGRkaGhobGxscHBwdHR0eHh4fHx8gICAhISEiIiIjIyMkJCQlJSUmJiYnJyco\nKCgpKSkqKiorKyssLCwtLS0uLi4vLy8wMDAxMTEyMjIzMzM0NDQ1NTU2NjY3Nzc4ODg5OTk6\nOjo7Ozs8PDw9PT0+Pj4/Pz9AQEBBQUFCQkJDQ0NERERFRUVGRkZHR0dISEhJSUlKSkpLS0tM\nTExNTU1OTk5PT09QUFBRUVFSUlJTU1NUVFRVVVVWVlZXV1dYWFhZWVlaWlpbW1tcXFxdXV1e\nXl5fX19gYGBhYWFiYmJjY2NkZGRlZWVmZmZnZ2doaGhpaWlqampra2tsbGxtbW1ubm5vb29w\ncHBxcXFycnJzc3N0dHR1dXV2dnZ3d3d4eHh5eXl6enp7e3t8fHx9fX1+fn5/f3+AgICBgYGC\ngoKDg4OEhISFhYWGhoaHh4eIiIiJiYmKioqLi4uMjIyNjY2Ojo6Pj4+QkJCRkZGSkpKTk5OU\nlJSVlZWWlpaXl5eYmJiZmZmampqbm5ucnJydnZ2enp6fn5+goKChoaGioqKjo6OkpKSlpaWm\npqanp6eoqKipqamqqqqrq6usrKytra2urq6vr6+wsLCxsbGysrKzs7O0tLS1tbW2tra3t7e4\nuLi5ubm6urq7u7u8vLy9vb2+vr6/v7/AwMDBwcHCwsLDw8PExMTFxcXGxsbHx8fIyMjJycnK\nysrLy8vMzMzNzc3Ozs7Pz8/Q0NDR0dHS0tLT09PU1NTV1dXW1tbX19fY2NjZ2dna2trb29vc\n3Nzd3d3e3t7f39/g4ODh4eHi4uLj4+Pk5OTl5eXm5ubn5+fo6Ojp6enq6urr6+vs7Ozt7e3u\n7u7v7+/w8PDx8fHy8vLz8/P09PT19fX29vb39/f4+Pj5+fn6+vr7+/v8/Pz9/f3+/v7////i\nsF19AAAACXBIWXMAABJ0AAASdAHeZh94AAAgAElEQVR4nO2dCVwU1R/A38zOXiy7XAoiIIcI\ngiKKKBrmjfD3yPv+561Z5vnPtPBIzTQtjw5PLMszrTSPtEjTrExN0zwqFVHzvhBFEIF9/5md\nBRGWZQd2Z4bH7/v5tDO7++a9X3zdOd6JMEA0SOoAAMcCggkHBBMOCCYcEEw4IJhwQDDhgGDC\nAcGEA4IJBwQTDggmHBBMOCCYcEAw4YBgwgHBhAOCCQcEEw4IJhwQTDggmHBAMOGAYMIBwYQD\nggkHBBMOCCYcEEw4IJhwQDDhgGDCAcGEA4IJBwQTDggmHBBMOCCYcEAw4YBgwgHBhAOCCQcE\nEw4IJhwQTDggmHBAMOGAYMIBwYQDggkHBBMOCCYcEEw4IJhwQDDhgGDCAcGEA4IJBwQTDggm\nHBBMOCCYcEAw4YBgwgHBhAOCCQcEEw4IJhwQTDggmHBAMOGAYMIBwYQDggkHBBMOCCYcEEw4\nIJhwQDDhgGDCAcGEA4IJBwQTDggmHBBMOCCYcEAw4YBgwgHBhAOCCQcEEw4IJhwQTDggmHBA\nMOGAYMIRQfDx3wE7cVz4X9/xgo8gwG4cEfznd7zgX1C2w8uoJGSjXwQfA4IrECCYcEgSfPeJ\n3QOp+BAj+PpQL6Ss/6kDgqnYkCL4vHf02hM/TnN62SHxVGBIEdy6ren8fFC5wwHhVGQIEZyC\n/uR3BnWxfzQVGkIEb9Obd1bUsnswFRsQTDiECIZTdEkQIhhuskqCFMHcY9Kf8JhUHFIEQ0VH\nCRAjGENVpUVIEgxYAAQTDggmHBBMOCCYcEAw4YBgwgHBhAOCCQcEEw4IJhwQTDggmHBAMOGA\nYMIBwYQDgglHPoLztm8qYBYIthfyEZzq6VaAE8pwSBmVEPkILswy9NDhZVQSxBeccSndWFoa\nEGw3RBVsPDq2pg4hpK05xvrcLyDYbogpOLsXQq7Rcd3iot0RejHHSkoQbDfEFDwdNTnAa809\nFIfesZISBNsNMQUH+GUV7OfUC7aSEgTbDTEFK3sWejNKZSUlCLYb4v6CHxfs5zYIspISBNsN\nMQXPeHoNPhyHZllJCYLthpiCn/Rm76Ibtese39gDob7WRoqBYLsh8nPw6EAN+xysCRx91Gpl\nBwi2G6LXZBnTL0JNlohAVSXhQFUl4VSCqsqs4ymlnjLIhfiqyrPtFQi5vJFVekoyIb2q8rRb\n/I9pFz/3aWvthEEypFdVtuiYx20uui21U4YVDcKrKi8h893c5Fj7ZFjhILyqcg+Tx+9s8rRP\nhhUOwqsqf6LNxazxsU+GFQ7CqyrvKXfzO4M62ifDCgfpVZUDw+9wm++ZXXbKsKJBerfZtKjq\nM7auGcYk2im/CgfpgnHW3CYGvw7f2iu7Cod8BF9LaFtAGHrgkDIqIfIR/Oi9uQV0hcYGeyGm\nYJdnsZISWpPshpiCV0YjFBBZgJWUINhuiHqKzolHW2xKCILthrjX4O0gWGzEFXxVt82mdCDY\nbsjnLrowINhugGDCAcGEA4IJRyrBV+E5WBykEpyKrOUCgu2GVIIzk5OtfAuC7QZcgwkHxiYR\nDoxNIpxKMDapKFnJi5cfrDSjlYgfm1SM7dXU9UMUkaccVoC8IH1sUjH2KhMzML7W1fOqo0qQ\nF6SPTSpG1EjTJif6FUeVIC8IH5tUjH/RSX5npZ+DSpAZhI9NKsYh9Ijf2aOoHPdZhI9NKsY/\n6BK/s9HDQSXIDMLHJhUPwWc+v9O9m4NKkBmkj00qxlInbjiacb7yiKNKkBcyqov+8/cC3nTg\nc/CbdOyrw8J0Gx1WgLyQj+DzNCqEA4euHJ/arf+7VxyXv7yQj2Ccca+ABVAXbS9kJLgQ0Nhg\nN0Aw4YBgwqnsgvf3jwh6YXWeSKVJQCUXPEPR64OkUfr2j0tPWkGp3IK/ZXZwmxSf10UpTgoq\nt+B2w/jtWj2xP+HKLdj1K357B1nvI1aBqdyCNebZs7LQb6KUJwGVW3C4uWnpMHVTlPIkoHIL\nnlnjHrcxdmojSnFSULkFZ0TW/T7jybEuLuT2sazcgvG9gQytQs+T67eyC8b44aE9xF5/OSq9\nYNIBwYQDggvI+XXFaut9ASsiIDifn2oqggNQA9Lut0Cwmd+1L9/F+Eq3qpfELtmxgGAzrXqZ\nNrlNB4tdsmMBwTxptPnvsJawEQ8gmOcMMj8N/4YyRS7asYBgnivoL37nWzVZN9Ig2ExN82jH\noYS1O4BgM6u0piW0PmH2iF2yYwHB+STSLSeObaTiVilN2fLFKVJO1CC4gKOT2nedehbj1FbI\ntSqKOCx+BI4ABBflpl+bUxhfHuB8QrIQ7AkILsroSFMPS2NXMu627CfYnhctKQVXT+K3v9D3\nJIvBjpRXcEFj+Z8t7BBNPhIKfkId4HfuIiLO0eUVXJsfSX1/nMKeJ24pf8GaHfw2BV2QLAY7\nUl7BmqBUjPNWe6IA21ZEssK1lALellBw6+H8dmF1Ip6Uyit4r873n6NNkHpauWtwzyORpnAo\nhd3MF9zmoGGhZCHYk3LfZP3q4kKhF1LsEMplefyC8ftMm6mzuzAjjfuHN2ky/CfpArEL5b+L\nPuqBlgs5Wv4Tgh8f26bZiD14oqLLnDldFBMljMQO2OEx6VS1ABt/wBVqQvDV2r3cZo/2M6kj\nKRflEZy/Tqw3citluVhzWfKYENxG6kzlt1PqShtHOSmP4NBnKfU4eUwIbiPpyFwZfahiLzdf\n6SYEt5Vr6G9+52/0Ws9eM8/dWTSo06R90sZUFson+NHRbX/n2nycLCYEt5Uc3VZ+ZzBq+PLI\nBoxzwKDxbRUDrF1YZEm5BCe5sndMUTYPjpfDhOC207+ZSeZmqhn7ellDb2c3Rz0nSRuUcMoj\nOBmhFv2CUXVbK+XlMCG47Vz07HAa49MGp4vsm/GNRjXlPtysSZc4LKGUR3BrtME0y7etVT5y\nmBBcAH83QwYDQpu4/ah399Hc6Sdb/b3EUQmlPIKr1eFez6BBth4ogwnBBXFh27afkekZP2jV\naXSL26m6SdqQBFMewagz95qDugs5WOoJwQWSpTQtovn8lG0a7m4yTfGzxAEJpVyCuxfe2BP5\nCMbtu3Cv8/zampYAmONl7coiR0BwKZzQvZTG3kW7KL5j/1iLmLVFv8/+8o3h809KEZltgODS\nOBCkrBPGBDakAqKcXD4p+u3pUENc//rUq7bXBohMuQR7d+YwbzrbMSo5CcY5+z5e8lMuPv7p\ngl2Fn5Eu7z70EN/36Xqf3d/nnihVdKVRLsHPUupxLs9iJaWsBFvkYH2koVQj3wrk6242q+Xa\nQ688gv94llKPWxmNUEBkAVZSyl7wL5pBf+Vl7Ah2MbcW5+i3ShtQiYjaLzonHtnWc0v2giOG\nmDbXFD3MH9RcJV0wVhG34/t2QgSfRqn8ToAvv83WbZcsGOuIK/iqbptN6WQu2LjKybwIQD9F\nhmn7uZNc66hh6IpgHozRI6QddJvbX6yMu85utuqtdV+QFBAslIdRtdYfZ6bVr8kN6ujcuaGm\nUccAZops+1CLL1j+vSqtkxhwB+Nhtc7VZ++zVisO5/3w7mtJFwu+Tdv/5Z+y6hQg8vKyFahX\nZUn4fcy+PGzu1kE9pwOz5NnvHo1SMVWQ31eSBGYZMQVb71WZnjipgHj5Cs7kZ//PWdEW1R16\n7NnvjPH+O7Pxran82Ah5IKZg670qb/XpWUBD2fZkvJNM/cDvpaBLOPPo3huFvvxCx3cQf9sz\nq/iREgG9KoVwpiVSInqgqVoyqer9UWqkQjFHC77uYZ4lL0PznSThWUJMwRWqV6Ul/nLrdCT7\nfeda9diH38s+bzSr+c39nBN9nQpm82g0z7wTImgwj0MR9xdckXpVWqBd+zw28t7Ohj47pldp\ns7DKNdOn/Rvmf9/SPBgCVyvWbCwZYgquWL0qi3ObPshtjJ8EKvQxi3Ji3+A/PouWdwqoGhjo\n7qNXerZe+seAcAN6eXm8j0/8WvZ5MHdpa0//jt88m1P2ghZVgromixG0mIIrWK/KYhxGfL0k\n/t50efFeb/6cYYYt8HbReih8UdREHR0/NyBcTw1ct260rk9uVpzbxC9WD1eNLpxReoxX4uZV\nA5jpIgQt8nNwxepVWYSTfMdKjLcauNdA82wtG9Fi3LbZw4+o3ng0ClFo67vHzHZ9Lob9Hzzp\ntvB1X1MlyM+aDYUyGhpqmtrkW0aEPrii12RVsF6Vhcl2+ZzfGWmaYamX+Z6xLp12Bp3Gse31\nj3FMNb2fcsmTwDFr6EN5KbvGB7hyhzz4Lfml58yZ3DuwN0W1k98f2EmEoKEu2nYm+5znNruV\npjaxn+mN3OYa3Q6vrY6xZiP6E89ye/06OvYSQlrk7s2+IHTZtDaTGtGmurtrXSlGidAZPr81\nPo6PGQQL4HF7w+hVH/dhzDfLCxVdF6+e4EFtxauCcB6zAx3B8/XT0lBjX3Qkrzql2IAPItTv\nUYM63z/KTaL0rOHbQU1+evxkCeWVaspgcxXHxyyV4KsVsstO3uqOgeF9fsx/+9uLkTUSloTN\nwz8p7+LQ0Yq7uH/wCz/S+guuG25Sg0fVMq53YRTDfO+wSReEdH8e49F1H7H7l1AUv35AYlPH\nhyyV4FSrnfTkKtgis31u5QaOwdM17fBx1RxFlOdoPDKiC3X/OnUg4uUuHh5cnewdv7f+oK5g\nz6STB1jdMeEqTvS/Hh86PjqpBGcmW3sKrFCCMxuFfPWVKiGaqfGKe8991SmPdx7/aKB87t1z\nDap5+4KOeif99pbaURk51IHbSINoVMubvTS7fnh9o38LEYZJwDW4/Dx4RcN3HKYQojsmMuxe\n3VB2H7Vgn4aG0OyeemQ6TkeH2qLpt7Mn0cyXO7jUuglirA4BDf724MnpI4/Sfk15fOp39tQ7\n/Plf7+JD9I5V9DX2qh0W2v/3U9wvdaNumSFiMr6mXTUoYqHvg/mK06KEBg3+ducEs5z9X419\nLqwfu3lT/wXzNffpBb//xY1bo93/cYDxIjLMxzhqpijRyKfBvzAVWjBepUyY98kwhWr86nea\naNq3aUi/sHDlKENCVuAneLwyLGLZEKpFLsaDB4gSjHwa/AtTsQXjE8Ojg9q/PzsuoJGXc++p\ngzz0YTU7r87DYR9j/F24c2g31S42UZ8RosQCDf6OpE001+Ejo0uQ6XbqRW583lf69H00++lj\n7xWihAAN/g7kIH3OtH1Y1TSw5RC9jhVbs3udvuzN1ygvcbrKQ4O/A3k3yrzTn5/GZLGi25L1\ng2nttC/mNXYRaRpbaPB3IBPr/G/4Iq59YlTM5CHv/sH+pPuFe7d5a0Ijl2p1pnAT6eVsHd6o\n0Yvz5wx9fcPj0jIrI9Dg7zj2GxSd+oUzb+Mzzqp2AxpSI8wPDr/5VenQP4J+05gSoUTuVRBS\ndUxwC7J5PjlhQIO/w0jRD1B+w95VaRZXo7nRw79Ufc30+VW3Idw9107nubVquuzFM/QewV53\nHvbxuuOQKKDB32GMaGZMNKzPw+9rGX56i2+Uph4hExrwQxOXaatq1+EH2nV76JrT8JPa0xwS\nBdRFO4ygZdg4U6ur64Jq8ZXOeW6buU299/jvH1JNNdl4F/ufT+dYjKfHOiQK+QhOUcljUQ67\n4cKdl+/s+uBraqT5k9rLuFef/D61TFMfjNey/zXsHsr+oy59wu2yIB/BxgPJBYwl4Rdc6yN+\newX9l9tkLB/CtP+KPTs3mMt/noaaqDNxsjrT6NW5OcaJzR0ShXwEF4aIU/SYhvzkWTPcvNnz\n0ekgr1imo1PzNPxmOP8A8b7eTb8CZxpWbFPWeAdnBjpmDDkIdhhXPfqmsRfe5crPQlteeRTQ\ndZNhJr4c0Qnf9u7K3jAbP1Mvi/TXbjIuVOsC/dNvJPg7pmYLBDuOY7V0TdtVc1qGL8Wo/NX+\nzJvs08Mp6g98KlzbOMFHswBfa05TGieElC2iNfXPOSYIEOxAnuyc88Za7tHIuK9BxCp+GoDw\nxRjn7p47+XPTwKYD09onTNrwyaR5P+Q5KAYQLA4dXjPvNJ8hbsEgWBxGdOO3Rt8kcQsGweKw\nQ3Me31/8YlwCc9n0PmNp5wD/Ron/mN78NaGur1doVP+FDpjwEgSLgzEhOKm635CuDNOWq8M5\nE+RM1QhVOSu5Z+VFCsagpxDlX8PrV7sXDIJF4mFvpKlv0My8ENoHsw+90ZrtGF+Pqst8h7cz\n2jEDap7foqzy0giP2/YuFwSLxbs11i3Yxt5SH6LO4aSqIaamhbN01+a4acPIVGo/xpMDFFdD\n7d7VEgSLRcex5h2/T/DAnuiEab/+COYRXfetdZ7s/hGk3zKprb2LBcFi8Xz+j7P+QtxlILpi\n2m8zDKWi4A+X1mb3LyCf1fMa2btYECwWfQfy2xzXzfjVBMV+0xv/4YZcp6hx25wfY5ysVP8w\noqu9iwXBYrFR/69pu0p3DyerYk3PxVuVMf/FvWtXvWRYgo2dwjwuuq62d7EgWCxuNgw5xE2B\nqH3v1CXc2Vs94m7eOn24+wV8ztUtIlE9b4CamR5aPxPfO36/5EzSjqcJLBYEi8OqQK4bg0uE\nVldHiZD7G0NpmqIpqh7X1e5YKIVohGiux4NCz75Efms5k6112C/rCFsdAgSLwhva2afu/dZL\nMXylW9vk2ymrvBMubh4/ZvERvnda3uHZI0a/Hewx78IMxt/t2O9jFZ9aymQJM/HY/WOvMUuF\nlAyCxeAovdu0nerVrKOp2SjFZWXRNB96XMYXVGtz2iZgvFh/s3gmV7T8MSu0VwQUDYLFYHwr\nfpulMz//4onFOug0norxO3Ux/o2+gfN8LUx2uTCY/70bgxcJKBoEi0FBY2GYwtxleEO1omlc\n2WvroIEYP6EOYNzpf8UzebmXeafXKwKKBsFi0CW/FiuENo9u+NyvaJoqmzAe3g/jR9yM4/EW\nlpIfk/+M3GWMgKJBsBjMiOB/uLeUCvN04oOKTXIXN5K9jfJ9gneoH+BHrhuLZ/KpJz9+6bHn\npwKKBsFicElr6uye06NOvwhTm28ys6tomi2qn/Bd98m3aw/FxhF+j4pnku41lvtnYhzrJaTT\nOAgWhQ3KLqt2f9ig6sm7kX6zd2wYoZxcPM041StfvKVU+3+1NNZgUcpe55bLdy9v4bxXSMkg\nWByO9Q5U1X31GsaZM6N11eJ3WEqzJc7LOSKqpip02AXLmZwdHKIKGSys+yUIlpK8FK5++sHp\nTIwzTqef4zpapp/Oyjidge8fPJbN2jldUDP58LSFs7YNgGDpuDVQh5Bb53AK0SG1uInTUJVu\nYabp1CglN61a7YYMQrVMi45vrccmii7LFPEgWDKuB0Z9felcX+Tz843JFKWPfm+IIooKSFT0\n7ETTNDX+m2hEL7pxbJb2dYw/ZF47fOPgK4oyLAUBgiXjv9Hsmfm86uPqb191WhbO3MB4KnrX\nTfUxPkmjncucLvjFaGs8YW+46UOpar4V8X2D8C5bIFgqMjTcvO+zIvF7IYtq3lK4f4pxovqD\n/6iNeJKXcoMx+CXNLT3DnZTjX50bzh+SW1341EsgWCpOI65B4b/D8H76pd4HURz74NQjeFRH\nBuMuvv5v4d6x9XGsJ9erNrHt4IHmYzpaqMEsBRAsFf+gq+zrwIH4B2ZUtyOo5RSM+wSM7ahk\nPfv5vo27N6+LY6pyI8YnJQzrZz4mwUINZimAYKnIdlnHvi4Myp3SIKnaXY3TZvZ8zayKUz7G\ns9zp7dnek5i/1RS3UFPMpI9q8DXYmW7rrWZpsRjZCL5Y3a0Ap0ogGE8IvI7xTcMo/Yr7VSbW\nVl7B+CVqhNYwAV9UoBUTq9yt5+fMTdeyXHX2jqtplQjjqOoZgkuRj+CcLZsKGF4ZBGfEVpv1\n3TfdkM/qfeNoplrVHlHqGBQ4SRMWrGEfevtP8EHM9B/X91WsxHi7pv3qfUktnfcLL0U+ggtT\nGU7R7N9+XkOtS+zMXv50za6dAyiVwtDinR5+tF5Pe3uqEeXSe1CIovoLpuFKJ/sE0EEDyjJG\nHARLSt7T17ynb7idvLyn3+Miu0IAwYQDgglHTMEuz2IlJQi2G2IKXhmNUEBkAVZSgmC7Ieop\nOicebbEpYaUQnLdzyotv/fD0/c0lI1q36TVnybQ2EYH1+oxPCK5Rr2VkVZ2Lb3CtqKZRIa0/\n+Lqdl3v4jPSr86KqGXwbzfhu5ehhC0roGvAUca/B20FwATeaatoMbqGMzx+JtF5XXadyUmiQ\nim8PpiiFebEt038eTuyLUkshlZIyfcAgbed+4cz8UooRV/BV3Tab0lUCwXkxMdwAhXPhHfj3\n+5mZnv0y8mpRAZTnpWMKyhmFN2MN16Eod3Yz0UDXp5D2Fr7qh5DWr92t+Qqk9hmC8UbV59bL\ngbtoifjGyTQTGv5LcdC0bTV4Wu0c/L0qjmLU6VMNlI9CoXDS0M4K7xQKKf+zjUbMc5pEjBtQ\nKNTzAcZBTsxI6jzGs/ytz74NgiVidH7H6MZvc6/Zih+em47xpDYDUEvtrqZuQShG76RjT9cK\nTS63+GGuGqGdVDROpxikHsoe4KZAsf4ruXkBrNdvwdqFEtF/uHnnhfHc6y10qvZSjIf374+6\n+a4N1UahHnoNg6qxF9vr3CX3oSuiUlEATmXf0FPYA9jzdmiTuRhnoYNWy4G1CyViQrx5p75p\n+ugc1a4Wb2I85flB6Dl1cnP3GqiB3smZ/bkqdHmsX6cn7G3Xl1RTnKFgb664RfE8GNSy+qdc\ns3Kq1XJg7UKJSFbxp9bDFD/esEO3uf6P8C9MM4p2zpyjowIpmlarFW601wn2GtxjjQIxdTWz\nMY5BKNzlJsYBGmYkw96mTQyzXg6sXSgV8bU5s7/6mdeoPKYZHdDuWlYoFYRcj+1naDWK6Mgg\nFEDTGhqhvmq6FYX09/BxD/YBKTD65Dj2LtptIs5ZwGy3XgysXSgV6d2osP8EU4PyV8RKrubs\nTClUKqTmZnLgnoMLFrDgdnXcUzH3aGww8F8wSNm4tYdhTSnFwNqF0vHH0skrzjx9m7ll1uBB\nr6/Zs2xoh9iEqUkjmjfuMLxzRPXApq1adBrYtfmAnSdfqRMcvz4nY1O3ugHRPZP+2Tln2hel\nTl8KaxcSDqxdSDiwdiHhwNqFhANrFxIOVFUSDlRVEo58qiorXcd3cZBPVWXhoSvOKKv44UBZ\nkGdV5S8ou4xlAEWQZ1UlCLYb8qyqBMF2Q55VlSDYbsizqhIE2w15VlWCYLshz6rKSiU4y8Ii\nHEZugpZHD+7m/xmu3sLZ926msU8hualnz2bfzMNpl/CtGzeuP0zP+/tffPPKVa5h+BI/L96D\n+4VmW5Jnt9nKIzh3QW0F8hv/7PLu+1vrkCHUh0KIrptkxMejFfzoBkRVNfB9PGgaFUXlwn6m\n6XM30ZOrK+x53pwXCJaUnM7u7/5yLCm0duEZzj5TDPr2WKQK+c6c4u3tNOw7Bd2liYpCVJvq\n3AyHGtPqLCpTPx46v1sPFczuaT/6qAVSudL9Ppzk5mb4nc9MKsFXrY4urDSCl7ie5TYPIgc+\n/eyK0wcYz/VqUaNqOr5bq7famfp2rdatnq8i7Lhp3Z1JOJD9PSuRriGlpbn1eBDdjfLQIUbb\nAeM4RHHL7lzxjgzLNeUmleBUZC2XSiO44VR+u0P9tPJ9Xih7jxI8hfrNYw3G61xa0f/Brdp5\nqHeEUAPZU3Uw8zhTyWqlKGV4axqpPShKpbyvoSgfHU1fy3V3ok1iF/kyP5lyk0pwZrK1qVMr\njWCNeeL3NPRHwWcDB3HjFRZocPxkjC+itvRs7PlCLEob7lyXPRcH0vg0ckGUQldF8T5Cdaso\nFAw605BCvb3dFXv+RZ78spdHUAi/vBJcgyXFkuBBnGCqkOB3sGfnAsFBTwUz8hRsS4N/pRFc\ncIrWPD1Fz3/2FN2aO0XHV9EUOUXTlk/RHuZT9GI/KU7RNjf4VxrBNt9kuVaEmyzbxyZVGsEl\nPCYNrpiPSbaPTao0gnHuQrtXdHhJVtEBDf4WIaeqEhr8JQAa/AkHGvwJR54N/keK3UIAZeaI\naIJtb/DHx39/iuZ/a2zDe4iNCZ+hW5iNCYd62ZjwdWVZ4lijG2tjwhov2pgwsh/757Ne42AR\nxzf4P4PO4kJ/FghdJjwc9uGtpY0JV1q77y/MbnVZ4sBuX9uYsN5iGxO2f71MgYhQF/0MILgI\nIFgQILgYINg6IFgYILgIIFgQILgYINg6IFgYILgIIFgQILgYINg6IFgYbt/bmDDik7JkP7ud\njQk/C7cx4V59WeLAXrb+Q45eamPCLlPKFIjYglNtXdTt3zI1QWXcsDHhk8s2JjSWuhSKRS7m\n2pjwyuPS05i49aBMgYgtGBAbEEw4IJhwQDDhgGDCAcGEA4IJBwQTDggmHBBMOCCYcEAw4YBg\nwgHBhAOCCUdkwStdbEn1aFI9p1qDrwnNvE+sieW2hfBkVpAqaGZJw+YKhWA9YYkhCDlMQAhC\nsxVXcE4jWwRnR6A6A55DLv8IyzxPzY/AS7QpBGNf5NvDB/WxPLaqUAjWE5YYgpDDBIQgOFsx\nBV/bmYBsEbwQDczF+DPUQlj2l9EEASEcRTFZOKsxOlZaCNYTlhiCkMMEhCA4WzEFc7Mu2SK4\nFbrObZ6jhPVS+RGV3r/paQij0QH29QAaV1oI1hOWGIKQwwSEIDhbMQVv27IlwBbB3gGmTR90\nQlD2SegHASEEuXITFOS4Wu5eWSgE6wlLDEHIYQJCEJytyDdZkbYI/sN08c3zotIE5f0GmhPl\nFDLkuk0hGDXRpnfRulJCKCVhSSEIOsz2EIRnK0fBJvLGoW7C8u6FqMZ9w5D7OVtCSEd8D9s4\nlGE9BFsSWghB0GElUzQE4dnKVfD1nsjnX2F5P6f/kv2LvIXibQnhIupuetcNXbIegg0JLYUg\n6LASKRaC8GzlKdj4sQE1Sy1TCbkhpayVmP8L5v8dxKF0y8nyQyg1oeUQynCYLSEIz1aWgu+0\nR55JtnYdL8qL6LANIRg1jcEDu98AAAJkSURBVE3vop0sP1EWhFBawhJCKMthNoQgPFs5Cs5s\ngjoKu7/ieHyd/+UORtYrSMwhBHpwgyxyPWqWFoL1hCWGIOAwy1gOQXC2chQ8FY2zdYRLIS7z\nlydjhNr6b98cwqumH/ohNKa0EKwnLDEEAYdZxnIIgrOVoeDc6m5luvNsRu9k/7jz0FibQjiK\n2uXinHaFZmcvIQSrCUsOQchhAkIQnK0MBV9ALjE8wpobTulQ6/4RKKKU2w9zCMbeKOrV+qh/\nqSFYTVhyCEIOExCC4GxlKHhvwcR9qcIyP9PLT9twamnLUeeHkD0jQBs713KzzDMhWEtoJQQh\nhwkIQWi20B5MOCCYcEAw4YBgwgHBhAOCCQcEEw4IJhwQTDggmHBAMOGAYMIBwYQDggkHBBMO\nCCYcEEw4IJhwQDDhgGDCAcGEA4IJBwQTDggmHBBMOCCYcEAw4YBgwgHBhAOCCQcEEw4IJhwQ\nTDggmCU0UuoIHAcIxiCYeEAw4YBgwjEJ/lETekPqQBwACMa84F91gVekjsMRgGBsEnzMxTdV\n6jAcAgjGnOBTHkLXeKkogGCWUO9qiDoqdRSOAQSzhCL/5ahpGWa4rQCAYJZQt/P4BbRa6jAc\nAgjG/F10itrzvtRxOAIQjM3PwdPKvoaVnAHB2Cz4kb/ipNSBOAAQjPOrKr9GLcq8gIZ8AcGE\nA4IJBwQTDggmHBBMOCCYcEAw4YBgwgHBhAOCCQcEEw4IJhwQTDggmHBAMOGAYMIBwYQDggkH\nBBMOCCYcEEw4IJhwQDDhgGDCAcGEA4IJBwQTDggmHBBMOCCYcEAw4fwfABcmrmG6/6YAAAAA\nSUVORK5CYII=",
      "text/plain": [
       "plot without title"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "suppressWarnings(\n",
    "    plot(degree.distribution(final_network), log=\"xy\", xlab=\"k\", ylab=\"Pk\")\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Use the `igraph::power.law.fit` function to estimate the scaling exponent *alpha* of the degree distribution:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "2.74535775518688"
      ],
      "text/latex": [
       "2.74535775518688"
      ],
      "text/markdown": [
       "2.74535775518688"
      ],
      "text/plain": [
       "[1] 2.745358"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "power.law.fit(degree(final_network))$alpha"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "## extra challenge:\n",
    "If you got this far, see if you can scatter plot the relationship between R (as the independent variable) and -log10(P) value (as the dependent variable).  When the effect size variable (e.g., *R*) can range from negative to positive, this plot is sometimes called a \"volcano plot\"."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAeAAAAHgCAMAAABKCk6nAAAC/VBMVEUAAAABAQECAgIDAwME\nBAQFBQUGBgYHBwcICAgJCQkKCgoLCwsMDAwNDQ0ODg4PDw8QEBARERESEhITExMUFBQVFRUW\nFhYXFxcYGBgZGRkaGhobGxscHBwdHR0eHh4fHx8gICAhISEiIiIjIyMkJCQlJSUmJiYnJyco\nKCgpKSkqKiorKyssLCwtLS0uLi4vLy8wMDAxMTEyMjIzMzM0NDQ1NTU2NjY3Nzc4ODg5OTk6\nOjo7Ozs8PDw9PT0+Pj4/Pz9AQEBBQUFCQkJDQ0NERERFRUVGRkZHR0dISEhJSUlKSkpLS0tM\nTExNTU1OTk5PT09QUFBRUVFSUlJTU1NUVFRVVVVWVlZXV1dYWFhZWVlaWlpbW1tcXFxdXV1e\nXl5fX19gYGBhYWFiYmJjY2NkZGRlZWVmZmZnZ2doaGhpaWlqampra2tsbGxtbW1ubm5vb29w\ncHBxcXFycnJzc3N0dHR1dXV2dnZ3d3d4eHh5eXl6enp7e3t8fHx9fX1+fn5/f3+AgICBgYGC\ngoKDg4OEhISFhYWGhoaHh4eIiIiJiYmKioqLi4uMjIyOjo6Pj4+QkJCRkZGSkpKTk5OUlJSV\nlZWWlpaXl5eYmJiZmZmampqbm5ucnJydnZ2enp6fn5+goKChoaGioqKjo6OkpKSlpaWmpqan\np6eoqKipqamqqqqrq6usrKytra2urq6vr6+wsLCxsbGysrKzs7O0tLS1tbW2tra3t7e4uLi5\nubm6urq7u7u8vLy9vb2+vr6/v7/AwMDBwcHCwsLDw8PExMTFxcXGxsbHx8fIyMjJycnKysrL\ny8vMzMzNzc3Ozs7Pz8/Q0NDR0dHS0tLT09PU1NTV1dXW1tbX19fY2NjZ2dna2trb29vc3Nzd\n3d3e3t7f39/g4ODh4eHi4uLj4+Pk5OTl5eXm5ubn5+fo6Ojp6enq6urr6+vs7Ozt7e3u7u7v\n7+/w8PDx8fHy8vLz8/P09PT19fX29vb39/f4+Pj5+fn6+vr7+/v8/Pz9/f3+/v7///95E5d8\nAAAACXBIWXMAABJ0AAASdAHeZh94AAAcrklEQVR4nO2dCZgV1ZWAT9XbeqdpaLYGbBahBRpR\noDEsUXZXUEigkd0BEyMiEhOFkCCiSJghmjghRNxnxMSYaAAzShLiOOoIIjIuAUFBkNWEHbqb\n7n71TS3v9Vtr30+f/9P3Xr+6devQf1fVrVO37gWOQA24HQBhLyQYOSQYOSQYOSQYOSQYOSQY\nOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQY\nOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQY\nOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQY\nOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQY\nOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOSQYOeYEn/vydNSiQAh7MCw4\nuv3ubvkAkNtt/odWBkRYi1HBdZMAigeMnjB6QAnA9HpLYyIsxKjgpXDVW5LWhvdGwwrrAiKs\nxajg8k41TZ/r+3ZXKvrh+4RFGDgZGhUc+nbSD3eGFUpuA8IytukWZXwPrm363HBFV4WSb0Od\nwW0QadTB27rXMSp4WeIcvHU0LFcoSYItw0nBFyfzreiBYyaOrWoFMOWiQkkSbBlOCuavg+/q\nksOfFXK63LVdMdlBgk1R/3Xis6OCBaKn96tnskiwCV65KgIlk/fGfnJcsKZUJQk2zkPBezZ/\n/NsRRe9LPzp8iNaYqiTBhtnOviK8Radd1iD+7KRg7alKEmyYeaOk92PB/xbfnRSsPVVJgg0z\nfGnsQ8Ua8c1JwcqpyvP/trKJuSTYKFdXP7R4/QnOFcHKqcrDY0c1cRmcM7iN5s7bhYEho0qL\nnnflEK09VfkrOGtwG82cPYUVTOebV68MbnCjkaU9VUmCDXJdqM832StKuo7Kd+MySXuqkgQb\n4wQ7soF75RsRCMBbsa+8maokwcZYDZuFt/oPAk03Cb2ZqiTBxpics176UAH7Yl85n6rkuXhU\nRTEJNsaNPa+XPnQuiX/lqOCL674z+ee1DQsiUFR9TKkgCTbG9O6BYb/n27HPwp3xr5wUfLa/\n0INkwk+h/cguUHZKoSQJNkJ0EctcHmIKr+8RDDfdMHRS8P1QvXX3Ssi9pZaLPgo/UChJgo3w\nSOHG6W2eWt0t74bQ401fOim4V2/h2rsKPuJfo/36KZQkwQY4m/8Md3FxHtsaCp5OfOuk4NxZ\nwutcEPNZ0/IUSpJgA/xXjvCLPfPuS3OHJn3rpOCuI4TXP0jn/1GdFEqSYAM83zH2YU1F0rdO\nCr4VnmqMf94WuEGhJAk2gLQH8ywZlvStk4K/agll08RPG2aGGaX+2CTYAGfznxXfa7slp/kd\nvQ4+NLdrL/HDVOj6ulJBEmyEHxds5F9PjC87mfSl05ks6W7SB58rp7JIsG7OLGwHACXjr86/\n7NPk711JVapCgvVyqm/3Jz/a8WDL9kteTe3sRoJxML/jpD5D7t11pOyBtAUkGAV1QSELzLC/\n+nnn9CUkGAODIVzQoqILMOvSe7ORYAwcAJj9m5d+2KILVMD51EUkGANVUPQ0/7anTZDtlraI\nBGMgBxaWH+XfHw/AyrRFJBgBFwDaBSMT/rZlITDpTwyQYAR8F6D/C98MMKHuUJa+jAT7n9/k\nFLPM3OMHO88ohD+lLyTB/mfY3d9vwzDQMgRwacZCEux/8jec7F3ejgUGSjK7MpJg3xMNv8Gd\nmJMDEJa6vadCgv1Pr1X8S/2uT/4SPJG5kAT7n4fLjgtv9cPHZVlIgv3PhapLf3f01Jbhbb7I\nspAEI+DMvFyAwLh92ZaRYJ+z/9G58359ov7v2y9kX06C/c2joYopE8tabpAtQIJ9zfpAt/xI\nvxWLI7LjjpFgP9NQEFny2uYV7YbeMFGuCAn2M8vgNeHt8CXjiuWKeHOsShKsjW4gPYT7bCE0\nyBTx5liVJFgTFxnoP37pAb4pDSVyZbw5ViUJ1kLdzZDT87bLC17hDsDNsoU8OVYlCdbCzDCE\nAIbOz/lgFPOJXCHvjFWZDAnWwEeQH4CCUChcEAzfJFvKO2NVJkOCNVAO1826pFW33NKS0KVf\ny5aisSr9yhsAwUiIzeUbrYX/lC9GY1X6lCh/OfKvT83LCc1goUqhHI1V6VN+CBCuWNu4u5ht\nrdhepbEqfUpInLau4N0nWDYkcyNJhMaq9CeDgAmz3QLAFIHiiGOu5KLPbN2tMnkwCVbhHYBW\nTHcIl0Ix5CnuL04KXvsj4fWzMfyhJbLgtFJJEqxCeygu+EF4Shehr+yziiWdFDxEWPFICfSY\nPrsSKmsVSpJgZb4C5jK24hcDw7zfUuWijgueBcsa+DPxKliqUJIEKzMJYGYhwwya1wFAZQJo\nxwV37y2eMqJ9BiqUJMGKPBkMwTW/n8Lw+y/0UCnruOC8W6UfpqaPVXnspuRpdc4Y3EZzoE64\n1AS26rOrh0Zgp1phpwX3u0r6YXjntKVnltzXxFjagxXoD93DHfmdN1wIMEqtsLOCF6758wrm\nZeHzRqhWKEmHaAVeAOjaA9iwsBv3ULnedFZwdUTMvhRy3NmJbMHfFUqSYHkaihj4hHumCPj9\nN/ekanFHEx2NB7asWzRpMH+tBMNk71ALkGB51jMd+97FcRceCoZzb1cv7k6vytqDystJsDxT\n81q8Glxey0VL2dAB9eIOCz6zM35QObxPoRgJlmVHCX911LGgcNDlDGjYgZ0VvOubAMwEae8d\npFQLCZZjS6RVJ7botpxJK5ezQaVkYBwnBR9qAYOr20LZfuEHEmyEmuJATohvP3dkLmXhl1rW\ncFLwLHiOb2gtgGHCwP4k2AjfgmfPld80i+UvRsJKcyIkcFJwd3E2kMZvwZMcCTbEXrbdY7/7\na4+Wo9jWzLjz6uU5ZwXnzRTfjhS2OUmCjRC9ht9xQ8yMdfPbdhigcR0nBVf2kR6geRzGNZJg\nA/BtVChlmZxruTtaK3bjSMJJwffBbGHETC56HdxzjgTr5jHh7hEz5lomsKmKVUoEJuOk4HOV\nAOW7+Q9fXwUtW5BgndQwwqDuEMwvZ7rCLK1rOXodXPez4e13CB8u/Lg9kGCdtAdmYlkBC3wT\nOtQ5fVBZWdx6ALzhi78qLCXBmTwD0I47XpUfEJ6//UzzavSEv19oxYvtePvB9XxLGn6ofTUS\n7BNWAhQEyvq23VXDn4jPqZePQ4L9wfssw3RkcnLbd5gMMF/HiiTYF5wtBAkhSdlCz5ok2Be0\nhHCeIDenF9/G0nAXOAEJ9gOTITR8dZ5wFcyqPIqUAQn2AccZ9tqx3PH+4kH6+/rWtUJwzaGa\nrOWMQ4JTmCjsu22e495iILdXo3r5ZEwKjm578Jpi/s+q+JoHlWb01gsJTuY2/hccAiZw9/4O\nwG7VubIpwQ3PD4BAv+o7F99Z3S8AA/9DbrQ13ZDgJJaIrWcmHBT249l61zYj+IMBBTPfiF9z\nn3tjVsHAHbrryg4JTvAVw5QVs7lTmZySucAqDXyRFTOC2/w0NaNybmUb3XVlhwQnKAUQJwfO\nEW4WPq97dTOCM5/hVnyqWwckuIkRAJFiCApPegP8RP/65hpZ57f/cZdlJ94kSHCcF4Wdt3h6\nJ0YYtfB7BiowJXid0IC+UuUJZCOQ4Bh7GWCYUAlU9iwqhpC2bnapmBG8GeDqW7tDhyzTMZmE\nBEtEi4ApAJZtC335HXiVkSrMCB4B68XRzX5mZMOKkGCJsUL/Z34fDgT5M7DS6K3ymBHcrrfw\n+qn2/kGaIcEi9wonYKaQ/5/3W2CssWNGMIwXXutBdkIIw5BggeelO4RMPpsbANDeSycFU4In\nJr9ZCQnmOS9msIIh6SbwfQZrIcGepYOQgs7pCMIsF3C50VpIsFcRbjEwLSBQXiSMDGj4fp0p\nwe3HC8TexhsNIQskmJsX66ITEE/DXxqux5TgVDSuTfMmaeFJ6XeawzeggwCPG6/IjOAdqaiv\nSPMmaeUIf1QugYDQvuIvkYxdAUvQvEleJJov7cAB6VX1iKcAzZvkRcr59nM/ho2dhv9hpioz\nghceT1tw7B7F9WjeJI08IrSbixfxl8BCFmu3qbrMCJ5beM+OxNGDP8EWKo/rQ/MmaWO3OB9D\nwSXL+wg3gU00sARMHaLfHAAV33n6nd2Hd7/z9O09oOot5fVo3iRN7G3NtAjzh+dISDhGX2Ky\nNnPn4Og7M0pjF0mlM/5XbT2aN0kL9T34a6PWUMyIv9ag2epMN7Iadz636v5Vz+3U0F2X5k3S\nwMkWwoWR8F+If2HTmzm6oXmTvEVdULo4as2Ko34rzFmntULTgrVkppKgeZMUaQ89YE9tCd94\nLugCsNp8haYEa85MJUOpSgXeBPiiw8+5i4PEBMcfLajRjGDtmakYlKpU4bVcvlnVqnAnxxWw\nkPOiFVWaEaw9MxXbFqUqlVkmXPdeWQCh2WtKdA3EoYAZwdozUxKUqlTmGbFzO8ytYG+4LABX\nW1OpGcHaM1MSlKpU5D+lC19gb2amLQPWoqdEzO3BWjNTEsp/EPvatmwirxnOm/Qy33LuGIg9\niATwpkXVmhGsPTMlofwH0fDqb5uY2/z24DMsFJVzNffGelI+ZVW9ZgRrz0xJUKpSnosDACLs\nNo5bxLAMlJlPcMQxeR2sMTMlQalKWb4qFnu3M08I43EAfGxdzaYzWZoyU/GylKrMzkHp3n4h\nMCt38afhGRZWTVO8e4ATjNR+hmLR82Ar66Z5k9znCCOMwQHicN+hAHSytHKaN8l1jjLxjsdi\nT448M13sMqF5k9ymhhU6x0aKY5K7ah7qWxtmBLdIRXU9mjcpCyfCUvtqULH4jNlcq+s3I/gJ\n/tqt/PImVNejeZMyqckT8xphyCsSro+KNdyT04epQ3T9WPiDjvVo3qQMjounXQak/hvAHLV8\nC+bOwRt0CaZ5k9L5tKl9xYjd3BUnVTaGOcGH8vV0OqB5k9LY2fTgniR6lw3boHmT3OP1VL2g\n2u/YCDRvkmv8T+z0G9f8N1u2QvMmucVP4ruvZJjZa89mzAoeFGfszPs/sCyqZiB4ZeLgLLwE\ntM5FqBezgof05qMTekr27R6Ga63Sgl/wA/HGs9TP3fwTDHKYFfzPK/ptPMOdf2PAiNoziyzq\nCNgMBI+It65iB2jjY3CoYVbwrHZS54MT7e7loiO1zlqsBnbBfSAVK6dDSMOs4I7VsQ/VPTlu\naZElMaEXfKN43g01XSVZn79KYFZwp5GxD8Pbctx3W1kSE3LBZ0skr/kxw+xBO7dmVnA1u0l8\n38RO4M53G2pRVJgFb286MEuCA7qnYdCFWcFfljK3PLr+sQlM8Z4jXeAli6JCLHhL2uk3aPn9\no1RMJzo+uUmMc+QObk+Hf7cqKryCX+d32TxWGAVa6uJeZvcGLchk7X993cY9/HuDdX1N0Ap+\nTUhq5Al6pT5YfW3fohWpylMfW9dPWwKr4CniqXcQFEmTqMAg+zdpWvCpJcIwLCWLTlkWEodV\ncGMPMXPF5JRDROhFCXMc2KhZwed6QrtbvjexPVQYmRJEDpSC343f3Q8GOuUK71Y9X6aIWcHf\nh/uFJ8pq79U5ra0yGAU/lGg5s6XdGYf8mhbcr6/UtGrsfYVFEQkgFDxMVNs3YfllZ7ZrVnDe\n9NiHafmWxCOBT3A4pjU3fpi27fZRGmYF9x4ovUcH9LEmIBFsguOdr/gL4EvED9fYMSFgVswK\nvgNWC8fo6Gq4w7KY0Am+OylzJSaiLcsIqWNW8IlOUDlv+bxK6GTlBHe4BOcl9DJigmONgxs3\nfR18aI7QIyEw55BlIXG4BNczkArjTPM5hgWZrLrdW3Zb/MQUIsFHIZ0vHN2+GcH1qVgYFR7B\nE9P1Bk6qr2Qlzk+rowU0gtMPzxY/3a0BM4KnpqJx7WY0GOmGjMPzVY7H4PB40c1qMNJhGX4f\ncz4ImjfJNjIOz2DplYZGrBB8eLy2OprVYKR/zNDLuBKHFYL3aOyL1ZwGI83J8NvBnUCcFNx8\n5k36R4Ze+I1LoTi7BzeTeZOmZejNdS0WJwU3l8FIAxl+1Ucgsg0rBDd8XZu1XDrNYzDSVzMP\nz5NdDIfmTbKYlpmt52NuxkODkVrKU5m7b4G7EbkxhMOZrbtV7kz4VXB+pl+tGVy7cFLw2h8J\nr5+N4f/ZkQWKc074U/CLmXqZGvXV7MVJwUOEFY+UQI/psyuhUqlh5kvBbKbfiNsxuSB4Fixr\n4M/Eq2CpQkkfCn49U6+uYQDtwnHB3XuLTaxon4Hpi7/8vImHfCe4OFOve8mNZBwXnHer9MPU\nvLSle1PuvvhLcEa/DR7l6e4dw3HB/WL3vId3Tl98+kQTq/0lOPPGIKjPA+cQzgpeuObPKxjx\nmY2NUK1Q0lfn4F5Zdt+R6qs5hJOCqyPiP76Q485OZAuUhnbzkeDPs+y+zB63o0rgaKKj8cCW\ndYsmDeavlWCY4tDI/hGc5drI+Y51SrgzGGmtyshBfhE8KYteO0b1NoFbo80q4w/BNVmOzhBy\nO6o0SLBhsh2dmc/cjiodEmyQ6mxH52+4HVUmJNgQx7MdnW0dc9IoTgrWPpGW1wVnOzpDN7ej\nyoqTgrVPpOVtwZ2z6XWn17M6jh6iNU+k5WXBK7IdneE+t8OSw9lzsNaJtLwr+GhWvazbYcnj\nrGCtE2l5VnDWky/zvtthKUCtaB1kdngWsHJ4IeshwZopzarXq42rOCRYI3OynnydHVDFCCRY\nE69l1Qvt3Y5LHRKsgcPZ9Xq47ZyABKuTXS/zqdtxaYIEq5H1ygjgQbfj0ggJVkZGrwd6tGuE\nBCsho9cXJ98YJFgeGb2MU0M9WwIJlkNO7wtuB6YPEpwdGb0OTHRkMSQ4G3J6PfO8gnZIcCbZ\nr3sBAi7GZBgSnMZBOb1+ajonQYJTWCOn1+s3jWQhwUlMkdVr3bybTkOCm8iT1as4nojHIcEx\nsvfWEPR661kjvZBgEbnrIr/cM5KHBHPcS3LHZmBUnoL0ASS4TF7vAceCsI/mLlj21AuMlRMi\nu0ezFlwru/P697o3nWYs+BJ5vT7NWmXDecEemTdJ/tjsz5yzHM1z3qQn5HdeCNq6ZcdpjvMm\nhRROvS1t3K4rOCnYE/MmHZfNafB6b7Nrq+7hpGAPzJsUVjg2M1ZOce0ZnBTs9rxJXyjsvJga\nzik4uwe7OW+SwpkXGB/2xdGIk4JdnDdpsYJdYOZYuzFP4aRg1+ZNUrjmxXtsjoF/3qQcxZ3X\nayMPWg7yeZN+qWQXmNHWbMXLoE5VKh6asR+bY+BNVSq1mvmd1xtTZtgP0lRlrqLdZrLzimBM\nVfZTtsv45+FeC/BOqvLic2ubmGpc8LNK6armtfOKeCdV+WXPrk20hQuGtrBVxW6zOfMm8Gaq\n8m2o01/9aRW7uG7ka8WbqUr9go+p2WWGGIza53gzValT8J9U7aJPWMnizVSlHsF3q9ltdu2q\nFLyZqtQsuFT5ikjYeSeaidb/eLPbrDbBEXW7zeqSNyu+FRxQt4usf6QxfCm4j+ppt5leE2XB\nf4KVbyJIsB/ZHqJf8Oa8SSmCw6EEGnZdYB8xGB9KvDlvUkJwSIPRZJgpBoPDijfnTWoSrGWP\nTbI7ymhkePHmvElxwWf12L3FcFiY8ea8SXHBQa122YcNx4Qcb7eiNTSYBbvnbA/Iv/hdMEPX\nu8p4W7DKIZoykep4W7BSI4sdbHsYGPC2YNlh1ynNrBWPC86W6GBLbd8+IrwumEtNVYbW2b5t\nZPhAMGEGEowcEowcbwrepprgIDSzTfev337B3IfvJ4jc+7xlVI6zrq4JFdbVdT9rXV1PwdNJ\nv733lZ/jzIoDgpPJ22RdXWMWW1fXA1dbV9dmCy/sa+BdkzWQYBESbBEkWB8k2CJIsEWQYH2Q\nYIsgwRZBgvVBgi2CBFsECdYHCbYIEmwRxZutq+vGpdbVtcLCgRDfzLOurouB7SZrcFjwvkbr\n6jpm4Si2549YV1f0C+vq4j43W4HDggmnIcHIIcHIIcHIIcHIIcHIIcHIIcHIIcHIIcHIIcHI\nIcHIIcHIIcHIIcHIsVfwxeVdw10fzD7gf/UQkbWG6zq3uE9en8Xn1bajuzLTgZ1eUFk4ZKmh\nwJTq0huXhK2Co1Og47fKoDo+KvzXTU/JreEaI9KnHxmsi6vrD5VTK6F/XZZlZiozHdjxchg0\nrSsMb9AfmFJdeuOKYavg7TCohqupgg9iP5+U/gaHdITfcwdgoam6uMfgjkau8TvwiyzLzFRm\nOrCZ8Bj/FzMFNugPTKkuvXHFsFXwXfAW//oWLEj9+lTnCVFuC78Xm6rr27CHf90Nk2W3Y6wy\ns4FdDFcKO+A/c27QH5hSXXrjimGr4K7FwtRK9cVpM99N63yC49bBn83VNQb28a/7YKzsdoxV\nZjawT2G2+H5Fkf7AlOrSG1cMOwVHcwaI7wPyU75+WYx0ETxyZV6P27R2dsusaxUs4l8Xwyq5\n7RirzHRgn8Ek8b0CzusNTKkuvXHFsVPwaRgjvo+G5BFGa7tcL7xNAqZqymVQssdoXY3fhREL\nhsOdjTLbMViZ6cAactsKn3aysFdvYEp16Y0rjp2C94M0z9EE+DLp20eZ/xPeBhf+jv+9PiAc\nFI3VFf21MCl46OmozHYMVmY+sCVw/d9P/6kLwF69gSnVpTeuOPbuwVIwo+F04suzraqTijT0\n0Dj5bGZdS+GWned23gzLs2/HaGXmA7vwbeFq5qbhcF5vYEp16Y0rjr3n4CrxfUBe0nXgWvhL\ncpnpsNVYXV+HKoR8QF3PyD+ybsdoZaYD47/668NLNjZUFWX/BRisS29ccWxtRXdpJTzI0NCq\nW+KraL8u4sMNtUekP8TZsNtYXW/DXPF9DryTbTuGKzMdWIyLJQNllxmpS3dcMWwVPE/8a3sP\n5ie+eg+kJ4oOSKebaGWkwVhdh0Bsq3HXwaFs2zFcmenAuFk3CpY2CUkTvYEp1KU7rhg2Z7LG\nNHD1Y2AHfzLZd0j86n7xUp5nKLuJD3cV3G2wrmgfZgP/9atMZcoy85WZDYybD2s57uilOSf0\nB6ZUl964Ytibi54MV87rB1P5j5tBmlvp8khs5vCP82HE1Eqo1NguyqxrRx4Mnf4NyP8wZZn5\nykwHdqwdjLy5GH7N6Q9MqS69ccWw925S3bLy3CErhfZLTPBhGBZf9umkTrn9f1xjvK4Dt/XM\n7fkvB1OXWVCZ6cD2T26bP3RT2jLzdemNS4LuByOHBCOHBCOHBCOHBCOHBCOHBCOHBCOHBCOH\nBCOHBCOHBCOHBCOHBCOHBCOHBCOHBCOHBCOHBCOHBCOHBCOHBCOHBCOHBCOHBCOHBCOHBCOH\nBCOHBCOHBCOHBCOHBCOHBCOHBHN7pGF62S4TP3Y7FBsgwbzgNhN5ru0EQc2juPgHEswLlsaX\na1igd5hAP0CCmwRzp6HU3UjsgAQnBEdzI5onBfANJDgheA8McjcSOyDBccHn37sq8IbbsVgP\nCY5fJgHkb3E7FBsgwfHLpHHdYdAZt2OxHhLcdA5uvB0edTkUGyDBiUbWTpjhaiC2QIITgo/C\nje5GYgckOFlwlbuR2AEJTgiugbZ17oZiAyQ4IZjrBStcjcQOSDBySDBySDBySDBySDBySDBy\nSDBySDBySDBySDBySDBySDBySDBySDBySDBySDBySDBySDBySDBySDBySDBySDBySDBySDBy\nSDBySDBySDBySDBySDBySDBySDBy/h8sJOXPwAA79gAAAABJRU5ErkJggg==",
      "text/plain": [
       "plot without title"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "inds_use = which(P_values_adj > 0)\n",
    "plot(cor_coeff_values_above_thresh[inds_use], -log10(P_values_adj[inds_use]), xlab=\"R\", ylab=\"-log10(P)\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## extra-extra challenge\n",
    "For each of the gene pairs for which R>0.75, see if you can compute the t-test P value for each correlation coefficient (don't bother adjusting for false discovery rate control).  Compare to the (un-adjusted) P values that you got using the Fisher transformation, using a scatter plot.  How do they compare?  Which test has better statistical power, for this case where M = 414? (If you are wondering, *general* advice is to use Fisher if M>=10; for very small numbers of samples, use the Student t test)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAeAAAAHgCAMAAABKCk6nAAADAFBMVEUAAAABAQECAgIDAwME\nBAQFBQUGBgYHBwcICAgJCQkKCgoLCwsMDAwNDQ0ODg4PDw8QEBARERESEhITExMUFBQVFRUW\nFhYXFxcYGBgZGRkaGhobGxscHBwdHR0eHh4fHx8gICAhISEiIiIjIyMkJCQlJSUmJiYnJyco\nKCgpKSkqKiorKyssLCwtLS0uLi4vLy8wMDAxMTEyMjIzMzM0NDQ1NTU2NjY3Nzc4ODg5OTk6\nOjo7Ozs8PDw9PT0+Pj4/Pz9AQEBBQUFCQkJDQ0NERERFRUVGRkZHR0dISEhJSUlKSkpLS0tM\nTExNTU1OTk5PT09QUFBRUVFSUlJTU1NUVFRVVVVWVlZXV1dYWFhZWVlaWlpbW1tcXFxdXV1e\nXl5fX19gYGBhYWFiYmJjY2NkZGRlZWVmZmZnZ2doaGhpaWlqampra2tsbGxtbW1ubm5vb29w\ncHBxcXFycnJzc3N0dHR1dXV2dnZ3d3d4eHh5eXl6enp7e3t8fHx9fX1+fn5/f3+AgICBgYGC\ngoKDg4OEhISFhYWGhoaHh4eIiIiJiYmKioqLi4uMjIyNjY2Ojo6Pj4+QkJCRkZGSkpKTk5OU\nlJSVlZWWlpaXl5eYmJiZmZmampqbm5ucnJydnZ2enp6fn5+goKChoaGioqKjo6OkpKSlpaWm\npqanp6eoqKipqamqqqqrq6usrKytra2urq6vr6+wsLCxsbGysrKzs7O0tLS1tbW2tra3t7e4\nuLi5ubm6urq7u7u8vLy9vb2+vr6/v7/AwMDBwcHCwsLDw8PExMTFxcXGxsbHx8fIyMjJycnK\nysrLy8vMzMzNzc3Ozs7Pz8/Q0NDR0dHS0tLT09PU1NTV1dXW1tbX19fY2NjZ2dna2trb29vc\n3Nzd3d3e3t7f39/g4ODh4eHi4uLj4+Pk5OTl5eXm5ubn5+fo6Ojp6enq6urr6+vs7Ozt7e3u\n7u7v7+/w8PDx8fHy8vLz8/P09PT19fX29vb39/f4+Pj5+fn6+vr7+/v8/Pz9/f3+/v7////i\nsF19AAAACXBIWXMAABJ0AAASdAHeZh94AAAgAElEQVR4nO2dCXwURbrAv56eIwfJ5ICEEI4Q\nwhEI4QrE5RA5QgQPINxyiQ9EDZeKIvoQBdf1gMV1XV0VH4vors9jcUV318XnqngBRkRhFQ9Q\nFOKBnBFIMjP9unqunp6k5+iePirf/+dvumdSXfXZf7q7qrq7CjiEakDvAJDEgoIpBwVTDgqm\nHBRMOSiYclAw5aBgykHBlIOCKQcFUw4KphwUTDkomHJQMOWgYMpBwZSDgikHBVMOCqYcFEw5\nKJhyUDDloGDKQcGUg4IpBwVTDgqmHBRMOSiYclAw5aBgykHBlIOCKQcFUw4KphwUTDkomHJQ\nMOWgYMpBwZSDgikHBVMOCqYcFEw5KJhyUDDloGDKQcGUg4IpBwVTDgqmHBRMOSiYclAw5aBg\nykHBlIOCKQcFUw4KphwUTDkomHJQMOWgYMpBwZSDgikHBVMOCqYcFEw5KJhyUDDloGDKQcGU\ng4IpBwVTDgqmHBRMOSiYclAw5aBgykHBlIOCKQcFUw4KphwUTDkomHJQMOWgYMpBwZSDgikH\nBVMOCqYcFEw5KJhyUDDloGDK0UDwRx8gKvFR7Hs/8YJ3A6Iau2Pe/YkX/A7UJ7yMFkI9vBPz\nNijYRKBgykHBFNL4U3AdBVPHixc4IGval75vKJg27rJev33fsyPTP/B+RcGUUWN5kSw8s4pd\nwncUTBmLRnuXP1jfEpYomCJ2rLv+kfLVvi89HhEWKJgavh/K9r+sOwxo8H5FwXTh3uAASJ74\nVZW9WviOp2jKmJOS9NrhV0ZlPsswBzmsZFHHVselV/AL98RBd8FofZtJdd+c8kRKg4JjZfxV\nY1eQ5Rewv7iDbh0dnpqlXVKBv1R0WSJ/lxIFx0r3P86aI6xkvdDlEb26KuunAmSUVVRVlGUB\nzG6USYmCo6dmaoGt+Npuf9icdZx8dd7DfC7+s5aCV8MFO7xaXTsr4G6ZlCg4av5im7Bp+x8G\nOsY0lFxYy3F7IG9eyN+1FFzQ4VxgvbG0SCYlCo6Wb1PuJwvXSHj5cFnyiBlpzIxzIQm0FGyb\nIvpSbZdJiYKjZW0vb431GGtd/Nf7xmVkbpMk0PYIPh9Yd/UrlEmJgqNl0iLfysC5IzLYbtf/\nLE2gpeA7g9fgXRWwViYlCo6W8df7Vob8mpfZRAItBTdM42vRA8dMqhyUDTCjQSYlCo6Wmy70\nLs87n2s6gcbt4MWdk/h2cFLnxTWynR0oOAJn7x2e020ar24v+7Lww+05Z5pOqXlPlufU19iT\npZRjfdutevbRqewGjluVtGbvsff+y/piM0mxq9J8nK6yM2kX8dXlZyzvc9zmIgBL+VvNJcau\nStPxIMMb5f+7ia9izSY//PxxXfOpsavSZLjnALRjLFXFDLuNe6hnxPTYVWkybrTAi3/NfKPT\n+HxLJbexS8T0BuqqPPJVgOdQcDPsZRk4/zHzzX7bEsjiFldG3MA4XZVfhrwS90ucZVDNRxOy\n+asv4+JK53CVU6HVl2mbI25joK7Kw8Ej+C5oplXXovm7Y/xzV/YCuJd7P2Xq5F9Bfrux7ogb\nGbOr8o8oOJyTbVZy3N0D01jrUe7DIfxZzrHyfOStjNlViYKbYFMOXzF5h13GMF1f+3mvDXLO\nRd7GqF2VKLgJbryUfI4cMA4YviUMHc9GtZUxuypRcBPccDn5/OmCVm1aJVuS/hDlVhoLPr33\nhG/t6CGZZCi4CTbmCRWYxmecxQueiOr0TNBU8GcXAjBV3wrr5XK5oOAm+Dnj18LywdQjMWyl\npeAjThg8PRfyvyZfUHB0eA5u/9TXq/u8dfarn782n/1TLNtrKfhKeJLj3MtgGGm9oeCoeK4A\n7JC1ztvefWdkEtiHvR5TBloKLhpKPt2T4QkOBUfF8csZ58Dlex5OX+L7wfWd3D2aptBScMpc\nYVGblnMCBUfDZ+0sF29e2y/z7Tcsu+LNQ0vBvUu8L7w9BJe7UXBkGnv1b83vMde1uSdH3hhv\nJloKXgHzvidLz1i4vg4FR8D95QNJN19E1s63e3jphHiz0VJwXW+AggMcaaxDphMFy9FwuxMA\nCkuEL7PmLZgRb0aatoPrN4zI20NWzq7KAxQsw8ftWbawtN+F8E/ybeG0zuvjzUmvF8BdB+Vq\n+y1csGuJBXJHLS+yvJ2VTrqshpRn/hRxo2Yw5hv+LVvwDyWQ3O/BKewKNvtV6Lj6qfmQuj3u\nzFCw4bjOAsBYSv/zN+sCuMI5dYiT7X8o/txQsMHwjAG2Z0VuekbbH+aN65lG+rHWR7z7JgMK\nNhb/nQYMdJo3e2y3rBufaXPxLObJfXJPRkQGBRuJutZAHmpn2FctN9m6bWuV26lKaZYo2EDs\ntEMKAwMYOxSsZ8ExK81S/KPSPFGwcdibDB25JBht72OBGx4GCwNTlD89jIKNwpH17VJhPtc3\n5Xcldge048/VKS+rkC0KNghPpXRlcsCybkNq6s2Z6WAFGH8i8laRQcGGoHYyk30RXMb2sT09\nIo2BFLAyj6uTMwo2AutYvu5shSw7u7D4/G3kDduy/SpljYL158zlfNu3z5//wqbarDZ4//hg\nBqSDIcUPCtadH4stjhdgin35KrBkWcmbd+z/qpc7CtadaYXWmVz6s6/bXu4OXcshmRmvrO8q\nFBSsM2d2soud1dwVozzzLluTVtY248aPVc0fBevK0z0ZgJ6tRnKfp897uMMKx8X22HXIg4J1\n5OxEtu+yddYxwLzH7exis9kh5+9ql4GC9WNXWxg4q9QCezpaMp/6pXF8CsyI7o3BWDCO4ENt\nMwOkwOmElGEoarOKh/OLf7A9PrMk25hkYKoTUIpxBLtefDbAAvqP4E/Gp1nSBh7l19Yywwsr\nh2Sy+fsSUY5xBIuh/xQ9ThhqJil7D8fVscP4mlbuWjUbR0FQsB4czwMmP7kqk8ktque4rD/Z\nNtcmqigUrAPVZAgG1sHeZ2GSXuROWUb3dCWsLBSsPWXAQI6TsdgdS8Ax6e3p1iz5kTsVgYI1\nZylA65QF7pUMpC61W1MYZsjhBJaGgjXmqa786dnmGMhxg8FusTAZjnUJLQ8Fa8sax2JILmoP\nlhe5pyGrG7APHU1sgShYS55sB2wmJCe/wFgcv9kCKXZ4INFFomDtcPcEvnplBci+ZR6wZEQ4\niPu97qhBwZrxXTZYk6HjAt6xbUFBOt9Ucu5NfKkoWCtOdoB0bsacZfZpAKk55PAdHnmsWOWg\nYG04uzITIOuxe7q5F5ZmCd2U1vs0KRgFa8Jf7cRpTtr4Vnd9DsOgqkfG81ocvhwK1obLARzk\noC1Pv8ZxiSUPijJrtCoaBSeeM2XAjPq0/WgGkscN2DuJr0ov/FazwlFwwllvEd4JLbg0mz9J\n2xsWMfCghqWj4EQzg2/7MlBYAcwNNu+EI7/WsnhjTm1Hj+AP0wWntmG5ZWxSa7AwTPewKX4T\nijGntqNG8CYGIKWM/+haPAP+1BUsv/qnxhEYc2o7SgTXDyZHb58+0IZl8zrCSiv003zCL2NO\nbUeH4FdYAJYB66gsaAUdyFMc/aOYB0dlDDS1nQgqBAuHL8MrbrWU4Rf8+iIdojDO1HZiKBB8\nNh/A2gmSSR06v+cIvp1UpF3jV4RSwT/5D8ozkSuHEaa2E2F+wZ9k8odvuyLo1Frod7YDjNEn\nEKWCYZNvZWXriNu1oKntZgqTODOdR8GVQM7OwD6hUySKBG/ZsgUWbBHY2Ds54nYtZmq7Woe3\nR6MHVDuFbg7I0e1VHEWCQyaEnRR5wxYytd0G7w7hLfdgLsoltedZ+gWjSPC2bdtg6TYvr0fZ\nwqN/arvb+BNyRU5rYpiFbGJ6p47RKL0GV/4rjjL37YkwM5uJBX+ZIbSOmPsDp7Y2p/SMR9Ou\nytrrZnNc3c18lZKdKfsyjnkFf8cyYIdUK2SP8fn9lb4BaSn4qzZwOeeZAHkz5veBjnL/rk0r\n+L1cKIE2zLhSK4zoJdxleEHniLQUPBk2urntcNkv/JX4AVgqk9Kkgs9Ns4CtKyRBaU4hCE9e\nzdC+b1KCloJzRvMfa+FTsu4p6yWT0pyCT1Zk3AQlb4KNsfTPh2QLgF6NXxFaCk4lszvdDt5z\n87Q0yV9rx44OUGzGIRyutJA+DUdtO2BZ4erLxFEBVR0tBQ91HuW4V7xTAZ3rOEzy17rVKwJU\nmu8I/lx4bpJvHHXZnQRMCnmsXffTM0HLvui/Q5e/1TeOK9zDcT9UwRqZlOY7RX/AQO5FeR1J\nr9WKw0Xk1tEdeofkRcu+aG5TKjgHDAco7GmDSrmmsPkEpwN79rFO3Crgz9M56VawHdQ7Ih9a\n9kVz3OnHhrXl94Ald8p22ee+TSb49DVd+LPy6a/Z7VwOMPz/oHW45k9uNIemfdECrqO1EUek\nMJfgfwmDXEGrvUtyXr8CLD3KBiVuyI2Y0bwvOipMJXg5QJIVUhm780y1JY2/+o76Qe+QROjR\nFx0ZEwk+wVec2awisDKZ8DvuQBaw6g3mrQb44Lsy/s2Q565sZER+FjKuLgZ4Tu+QQlEs+Pnp\no3yoFpOJBB9mgcnr8ck4cpfbQe78Mrv1DkmCUsEbAVKzvUTczhmKTEqTCG6808l77ZbBNQ6z\nZwt9Vz3r9I5JilLBPVPfiHpuzMfLAAr6BJBJaQ7Be5KF5oMFnuH+xtogF/5ooNqzH6WCHdfE\nsGFjJWyNKqEpBF/tbR7y7V7Lb/eDxcL01TuiplAquP11sWy5jR7Bro7CwUvetHK0Ij2U0NcQ\nfc9SlApenR/L7PJHUl+KKp3xBe8jN4zI66DJvOWUqWnAaPpSaPQoFdx4ZfenPj9+gqBeUMYX\nfKv39MzafadpsKg9mYZaKBXsdAb6KtULyvCCxxKnDr5aZfM/WKfmP29VUSp4YRD1gjK44PNt\nvadna3Iek2pxAHQ1VudVCNiTFTO7hfcWbJAEFhZSGQbaJ2YwfnVQQXD9/nd/jLotHB1GFjzS\ne1Lmq1WsUImG1sf1DkkOxYKPzk0C2PpyxSeqhcQZWXBDK3+VI1eoXAH01zskeZQK/qEISq6A\nre85Mr9ULyjjCv4n46s1k6VQhZZ7SdIIKBW8BNa4D8FWroa9Sr2gDCv4HpBgUXcmyQSgVHCn\nfh6OCOYGdlYvKKMK7iRITWKsfr/ZBux8lqBUcMpszit4Top6QRlT8AnWf9iCz/BKvUOKAqWC\nBxa7BMGesgHqBWVIwW8w3qsvr1m4NQgOo3ZehaBU8J2w6BwR/Bjcol5QRhQ8w3/0kn5JstZH\n5aZhglAquOFX0PZiGD0IStScGdVwgj1t/JddbwUaYKjeIUWJ4nbwufvzSW3jNlVfJjKa4B+S\n/cev3++teocULWp0VZ7ep/b4mgYTXENavEnCI+2+w/grvUOKGuyLjky+XyvjW2lr/NZRACWC\nc0NRMSojCT7vbx0Frr/z9Q4pFpQIHkLowP+D7t8OoFLN+ckNJPiA97xsC1iGx/QOKSaUnqI/\nSh9OZnfaP7K9mpcl4wjewx+/fcuFZ9t9t/bVbC1ogFLBk/O81ecz7aaqFBHBMILvJOfk1qyl\nXLhzxB/Mw/WOKFaUCg54nZqnMJJTtxnuDf/63t5L7yxgvN3QsELvkGJGqeC8wb6VCzoqjOTH\n6VMCDDDEGB17fKdlhr3EKlyJW32md0ixo1TweHhWWD4T9fvB0WCIU/SYQLWK7Q1L+EWFRpOV\nqYpSwftSoOqxlx6dCGkH1AvKEIK7E7XFwRYSlOgdUVwo7uh4o4+wA8pVvbWiv+BjvtZveeAw\nnqtzRHGivCfL/f6W+5/5UN1bK7oL/nPAa2fvwqHLePwqgF2VTTEseGLOFzqvppnx8iugVPDM\nIOoFpbfgYLeVb8i6J/WMRhmKx8ny07GLekHpK/gESGDf1C8YxSh++YzQcGRr6ahf1AtKV8H/\nkvq16Tqgt1LUugYf73CD8mAC6Cj4Ialfu2kvvwKqVbKuUdpVKUY/wRlSv9IxU82GaoLnJSmO\nJYhegt2M1K/ckKmmQCXBrleTSlWIxo9Ogt8QmfXeBjb8iwsRUSo41QvfsHhKvaB0ErxAevgy\nhhx1IzaUCr7Ex9yX1YtJJ8GpUr/SMelNCfZk+bFI/ap5zdEPLUd8jx7tBddI9cJmjSNIEJqO\n+B41mgu+Rao3RW4GejOh7Yjv0aK14KFSv3IDaZoL7Ud8r/vGcJNTJkn9FmhZemLRdMR3T83S\nLqSumtxlyUeyCbUVHNa7Yb5H65pHyxHf66cCZJRVVFWUZQHMlrvIaSn4i7Dq1deala0BqjST\n9m9+LZqnIFcHp3jfWQF3y6TUUHCVVC+jVcnaoEiwa92oh/nFzfxuaRPFPdOCDsG5khpLi2RS\naic4R+o3Q6OCtUKJYFcFMA9x3POQees19uTIU43Ypoi+VNtlUmomOOzyW6FNudqhRPBTMIGM\nwTkcXiWWI49IUtAh2Lfr6lcok1IjwWfD/N6mRbGaokTwCCfxe9JCRqtwt4s80/WdwWvwrgrZ\nEcS0Efxs2OU3wszzZkSJ4I4Xk8+/wQayGBb5hn/DNP4SN3DMpMpB2QAz5Ebw1ETwrVK/bOLL\n1B4lgh1zyOcSIO+PciOi6Mny1CzuTDoVkjovrpHt7NBCsFPqt3vCi9QDJYK7juA/3EXZwkNL\nnYqj29hz6mtD9GSFXX7vTXSJ+qBEcBX7KZnweTJZfw8mRrPp6b3+odGPHpJJlnDBDWF+dyW2\nQN1QIngn5D38l67wAr96oCe8EXnDzy7kKzJV3pdAyuXuKida8Fth1Ssjj+mtCEUdHevJg/8X\nu7izpdZoRiY54oTB03MhX+gK1FPw2DC/iSxNX5R1VX6y4YYtfMPnTNrQJ6J4++xKeJK/Zi+D\nYeSiraNgu9SvI4GF6Y0qfdGe6J4NLxKG/3NPhie4JgQ3bHk0wMxECg67/Jrzxd8o0fKZrJS5\nwqI2LedEE4K/6VYYoHUCBYf53ZGwooyAloJ7l3hHiHsILnfrdYreE3b5NfebKRHRUvAKmPc9\nWXrGwvV1+ghe3oKqV160FFzXG6CADOXx0wWQ6dRDcJrUryUhxRgJTZ+Lrt8wIm8PWTm7Kk92\nKrwECQ579JnO3skQ9Hrw3XXwdZm/JkZwWPXqdwkoxGi0oDcbwvweUb8M49FiBO9ocdUrLy1F\ncEvqnQyhhQi2Sv1SeXO/KbQU7AxFJqXagsOqz1T3ToagpeDHy/h2cJ8AMilVFhxWvXpCzdyN\njaan6MZKMgleFKgrWOqXieItG2rQ9hq8TQfBX7XU6pUXbQUfSX0pqnQqCp4u9Ut/72QItNei\nw56dtKmUsVmgXDAr9avmLMemgG7BYc2jltD7HArVgsOqz8fVyNVc0Cw4zK8KeZoOegXvbNnN\nIz/UCp7QwptHfmgV7GjpzSM/lAoOax4NVCcw80Gn4LDmkaqzOpkKKgVj9TkIjYLRrwj6BDdK\n/bbQ6rMP6gTfK738tpiHc5qGNsHZ2DwKhTLBYa/+UjRubHzQJTis+XuXyoGZD6oES5u/jN6z\nEBsAmgRj86gJDCT4kw8C3BqPYPTbFMYR/GWIoGjGnw4Fm79NYhzB3JnjAX4b8xF8EJu/TWMg\nwSJivgaHPRzb0pu/AegQTOesdKpAheCwu/s3Jiow80GD4LB3Q2sTFpj5oECwtPsKm0dizC84\nrPsqgXGZENMLxuavPGYXjH4jYHLBUr/YvSHF3IKlfrF7IwwzC14jbR6pOYUxLZhYcHup35yE\nx2VCzCu4ldTvuoSHZUZMK1jaPYnN36Yxq2Ab+o0OkwrG7sloMadg9Bs1phQs7X7G7qvmMaNg\nqV/svpLBhIKl3VfWhIdjZswnGLsnY8J0gqV+kxIejLkxm2Cp3/SEx2JytBdc942CGcClfjsp\nCqUloKlgT83SLuQJ1+QuSz6STdiM4NslemFVvJG0HLQUXD8VIKOsoqqiLAtgdqNMyqYFl0q7\nN1rSyO3xoqXg1XDBDq9W184KuFsmZZOC20j9xhlGy0JLwQUdzgXWG0uLZFI2JTgd/caDloJt\nU0Rfqu0yKZsQnIJ+40LbI/h8YN3Vr1AmZbhgvP0bJ1oKvjN4Dd5VAWtlUoYJlg6ugrcXokVL\nwQ3T+Fr0wDGTKgdlA8xokEkpFSy9vY9+o0bjdvDizkm8n6TOi2tkOzskgqVvl6Hf6NG8J8tz\n6uume7K+CjkN14n+Ir29j7cHY0CXvuhj7zYxN7Nnx/YAD4CoD0PqF28PxoKmgusfnFp5v9uz\nhj9We++VS/iOSLD09j7eHowJLQXX9SeCFv0P9Jg/lkn9TialSLDUr1zzGQlHS8G3wDUHDq0E\n20S+Av0izJdJGRQs9Yu3f2NES8E9S938pXY4fEK+lBfLpAwIlvpNjbPolouWgpNnkc+l8AtZ\nzJJ7U8wveK/Er9ys4UiTaHoE9+GPYO4i7xE8pKdMSr9gyQHcJs6CWzLaXoOv/eLQLcBW8dfg\nl2CBTEq/4NAHODrFWW6LRtNa9ACi6dpHoceCcdHVokMEXxhnsS0bbdvBv5syhm8H324DKJF9\nZqepU/Tt8ZbastGlJ+vHd47IP3fnF/xT0O/bSgttoRjzsdmwZhLe/o0Xgwv2dUSj37gxumCe\n8XLVMSQCJhCMKAEFUw4KphxjCt4NiGrsjnn3J14w99EHsmyBx7co53bYpEIutzIqZLLlJrsa\nuSxLDd9X8m+BNYkGgiPwQRzT7oSjzoXgdVXacH9XZRjGrRlq5IKCQ0DBCQAFNw0KDgEFNwMK\nFoOCEwAKbhoUHAIKbgYULAYFJwAU3DQoOAQU3Az6C/7YclaFXHZb5Ub6iZa3HSpkwr2myiPf\nr6gzbYX+grmvjJOL56AaubgPqZGLS5VcjCAYSSQomHJQMOWgYMpBwZSDgikHBVMOCqYcFEw5\nKJhyUDDloGDKQcGUg4IpBwVTjl6CH/fdFG9YW2gvXNMgWY01l+lDBB6NI5dfVpSmdJ13VGEs\nolwUxPLd7KKU3ivPKIwlFJ0ENw70qvHMgPaT82G6J2Q11lzcvlkgbos9l/re0GvOYHAeUBSL\nKBcFsRzNhIvmFsOARqX7RYwugo++crFvJMMaKD/HnRsEH4asxprLYbjB/2OsuWyAuS6O2wzD\nFcUiykVBLFfDExznmgYble2XUHQRTKbE86pZDDv4zx2wLGQ11lz+DY/4f4w1lxFQSxaDmdNK\nYhHloiCWwnwyPOROWKhsv4Sii+CXtm4t8KopzCDPyjVmFIWsxprLRnjN/2OsueQVCIvpsFdJ\nLKJc4o+lsacwwOsBmKZsv4SiVyWrj6DGk1QmfCtLFa/Gmgu3En7TP6XbVbVc7LnsOUA+3bnM\nCSWxBHNREouXe+FBxftFhL6CT8EY4VsF1IlWY82FmwrMoBnFkPUFF1cunHsZVCmOxZuLwli2\nLrwAJp5XHksQfQV/DZOEb1XwjWg11ly4wWnP8/v3Dqjk4sqldgrkf6s4Fm8uCmOpBki+z6U4\nFhF6H8GVwrcKOCVajTUXH65ucCaOXDx/SIehh5TG4s9FWSwcd37vBLhe8X4Rofc1eJDwrSzF\nI1qNNRc/s2FX7LkcGwc5G11KYwnkoigWgXN5jgal+0WEvoK5ztmkZeDK7hKyGmMu52u9c6vN\ngwMx53L2Arj0BKc0lmAuCmL5cOY2YTkKvle6X0ToLHgR7OJI029JyGqMuRz2XqY8vR2umHNZ\nBcvcvlUFsQRzURDLlzBX2LSz06N0v4jQWXANjHFxjWNgT8hqrLkMtbzC75j7YGnMubjaZQYq\np/HHIs4l/lg8hfYP+M8HSDtY4X4RobNgzzTov6gvzAxdjTWXfakwcmZv6H0q5lwOgrPcy1EF\nsYhziT8W7lXGOmZ2P2j3veL9IkJnwVz9nQXJQ+5pkKzGmst/pnZIHrDqXOy5vB4YJPCQglhC\ncok7Fo7bNbZ9Sp/lJyWbxrFfROD9YMpBwZSDgikHBVMOCqYcFEw5KJhyUDDloGDKQcGUg4Ip\nBwVTDgqmHBRMOSiYclAw5aBgykHBlIOCKQcFUw4KphwUTDkomHJQMOWgYMpBwZSDgikHBVMO\nCqYcFEw5KJhyUDDloGDKMZHg3YH36Ku5AfAP0V8Wwgl1i9rYIfmLuDYMjcsImEpw/iQvj8cu\neBtsib6kWrbNkh9jjM5bAApWwG4y/IyPI4fE88KrLfhtuDuWwEQFhMZlBEwqOBS1Be+A30ef\nOJ4CNMSkgolS96ZBzqwL/yl8+3lt/5ReG8lfGtaWp3a+npxgF2a7l7Xymaokl+6fAj/tmdze\nnj+xhiRyNt7RMankCX41kN9MkvhTjjt+bUlqv+W/BLOqdp5f1r31hO9/ubao1YhPSL6BjPwF\nkH9qou2CueuFiQWvAef4KSmWN8m3Ce2rr06Bv3Lc+cHQY1Zf6FpLrNwBmU970/9rKSzYdM7/\n0xdOduycEnB+SxTMy79uYSq8wAXze3MlzNx0kjvSEcpm94YeJwNZVaeO7XfThdBnYM/lFdDV\nxYky8hfAxyXeLpi7XphKsK+S9XtBsCe70xmOexOuJN96HOO47TCL49ZBtYvz3Cn8asl/K7Cx\n9wzq+2kVPM9/rofNZNPu/NH+BkznRPl5T9HzYT3HeW6GVYHtquGSRs4zEIae4zyj4WBIRr4C\neMHi7QK564apBPuYK+zIekthI39efW8/+UYOVHdSJcfltyVDkLl7JTfwvz4W3Ni//4WfXnuc\njJL/KmwgvzzFr3pSR3Gi/ATB9WwJGX/yXNs2ge2q4V3+80Z4hf9cC7tDMgoIDt3On7tumEpw\n6Cn6Mui1YZ/b++1zsnBWcqfh4kOEK2A//+tnwY39+9//09ldD5R4BQubZvMKgvkJgg/AIiHl\nRDjp364afuA/bwMygP/9RLAoo4Dg0O0CueuFiQWfuaUtQNulx8i34+RHXvC+QGfIu/yvx4Ib\n+/e/8NPJG3qxltKxXsHCpkRBMD9B8L9hrbDldbDPv101X40igkkniCBYlFFAcOh2gdz1wsSC\n+fPpB+v7Qz+3v5nEC/4ZRoBQHCYAAAGWSURBVG/18mNo4yl4ieQZDwteruPe8woWfvEq8OcX\ncgRPgp/9qcIEizIKO4JF26HgqJAI/mr1//ErnpF8ZScomMsqF/7+/sue5gWfsQtjOv8lVLAo\nP981uDcZY/18uyyuOcHijETX4LDtUHBUSAQfgrJ6fj8OYM+IBd9G5g3jahyjOangjZz/p+Mw\njFdwuDvprwoqEOXnr0Xz3tzL4ZZmBYszChYQvh0KjgqJYM8l0O2qy7LISOgiwad7waC5g9iM\njyWCXwcy66Pvp9FQOL3Sdqm1zXqRAlF+XsFHOsCg2SXe9mwzp2hRRsECwrdDwVEhvQafXNkt\nOav8cVeIYO7szX2TC678gpMIrq9Kyg5cE3+an58+YpNnfc5NYgXB/Hxdlcev7ZXSZ3kd17xg\nUUaiAsK2Q8FIwkDBlIOCKQcFUw4KphwUTDkomHJQMOWgYMpBwZSDgikHBVMOCqYcFEw5KJhy\nUDDloGDKQcGUg4IpBwVTDgqmHBRMOSiYclAw5aBgykHBlIOCKQcFUw4KphwUTDkomHL+H75y\nWks5sULnAAAAAElFTkSuQmCC",
      "text/plain": [
       "plot without title"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "ts = cor_coeff_values_above_thresh * sqrt(M-2) / sqrt((1 - cor_coeff_values_above_thresh^2))\n",
    "P_values_studentT = 2*pt(-ts, M-2)\n",
    "inds_use = which(P_values > 0 & P_values_studentT > 0)\n",
    "plot(-log10(P_values[inds_use]), -log10(P_values_studentT[inds_use]), xlab=\"Fisher transformation\", ylab=\"Student t\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "R",
   "language": "R",
   "name": "ir"
  },
  "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": 1
}