scripts/beastAnalysis/eigen.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Human                  0.001054\n",
      "Chimp_pts              0.006239\n",
      "Chimp_ptt              0.006782\n",
      "Colobus_Monkey         0.009223\n",
      "Mandrill               0.014447\n",
      "Sooty_Mangabey         0.015430\n",
      "DeBrazzas_Monkey       0.017116\n",
      "Western_red_colobus    0.017693\n",
      "Gorilla                0.018284\n",
      "Red-capped_Mangabey    0.066527\n",
      "Red-eared_Guenon       0.070091\n",
      "Vervet                 0.091364\n",
      "Tantalus_Monkey        0.092028\n",
      "Patas_monkey           0.168551\n",
      "Baboon                 0.188008\n",
      "Sabaeus_Monkey         0.217163\n",
      "dtype: float64\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Library/Python/2.7/site-packages/ipykernel/__main__.py:28: ComplexWarning: Casting complex values to real discards the imaginary part\n",
      "/Library/Python/2.7/site-packages/ipykernel/__main__.py:30: ComplexWarning: Casting complex values to real discards the imaginary part\n"
     ]
    }
   ],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "import pickle\n",
    "\n",
    "actualRates = pd.read_csv('../beast/supplement/discreteTraits/results/discreteTraits_2016-10-26_actualRates.csv', index_col=0, header=0).T\n",
    "\n",
    "s = actualRates.sum(axis=0) # Series of the sum of each row\n",
    "maxRate = s.max() # Max row sum value\n",
    "\n",
    "def normalize(x, total):\n",
    "    if total == 0:\n",
    "        return x\n",
    "    else:\t\t\t\t\t\t\t\t\t\n",
    "        return np.float64(x) / np.float64(total)\n",
    "\n",
    "for host, row in actualRates.iterrows(): # Iterate over rows   \n",
    "    total = row.sum() \n",
    "    extra = maxRate - total\n",
    "    # Add to diagonals so all row totals are equal\n",
    "    row = row.fillna(extra) \n",
    "    \n",
    "    # Normalize so each row sums to 1\n",
    "    actualRates.loc[host,:] = row.apply(normalize, args=(maxRate,)) \n",
    "\n",
    "eigenValues, eigenVectors  = np.linalg.eig(actualRates.T) # NOT SORTED\n",
    "rearrangedEvalsVecs = sorted(zip(eigenValues,eigenVectors.T), key=lambda x: x[0].real, reverse=True) #thanks stack overflow\n",
    "principleEigenVector = rearrangedEvalsVecs[0][1]\n",
    "total = float(principleEigenVector.sum())\n",
    "\n",
    "principleEigenVector = pd.Series([ float(x) / total for x in principleEigenVector ], index = actualRates.index.values)\n",
    "principleEigenVector.sort_values(inplace=True)\n",
    "print principleEigenVector\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "principleEigenVector = { host: score for host, score in principleEigenVector.iteritems() }\n",
    "pickle.dump(principleEigenVector, open('../figures/sizes.p', 'wb'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}