eigen_byarm_tutorial.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "name": "eigen_byarm_tutorial.ipynb",
      "version": "0.3.2",
      "provenance": [],
      "collapsed_sections": [],
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/gspracklin/genomics_scripts/blob/master/eigen_byarm_tutorial.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "metadata": {
        "id": "qt4vY1jT5z_-",
        "colab_type": "text"
      },
      "cell_type": "markdown",
      "source": [
        "#Eigenvector Analysis by chrom arm"
      ]
    },
    {
      "metadata": {
        "id": "BqnJsU0sDddL",
        "colab_type": "text"
      },
      "cell_type": "markdown",
      "source": [
        "###Install requirements"
      ]
    },
    {
      "metadata": {
        "id": "Ld3QU93H5sE3",
        "colab_type": "code",
        "outputId": "c632680a-6004-401a-f344-fbff7544afe2",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 1516
        }
      },
      "cell_type": "code",
      "source": [
        "!pip install cooler\n",
        "!pip install bioframe\n",
        "!pip install cython\n",
        "!pip install git+https://github.com/mirnylab/cooltools"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "Requirement already satisfied: cooler in /usr/local/lib/python3.6/dist-packages (0.7.11)\n",
            "Requirement already satisfied: pysam>0.8 in /usr/local/lib/python3.6/dist-packages (from cooler) (0.15.1)\n",
            "Requirement already satisfied: pandas>=0.19 in /usr/local/lib/python3.6/dist-packages (from cooler) (0.22.0)\n",
            "Requirement already satisfied: pyfaidx in /usr/local/lib/python3.6/dist-packages (from cooler) (0.5.5.2)\n",
            "Requirement already satisfied: multiprocess in /usr/local/lib/python3.6/dist-packages (from cooler) (0.70.6.1)\n",
            "Requirement already satisfied: click>=6.6 in /usr/local/lib/python3.6/dist-packages (from cooler) (7.0)\n",
            "Requirement already satisfied: cytoolz in /usr/local/lib/python3.6/dist-packages (from cooler) (0.9.0.1)\n",
            "Requirement already satisfied: numpy>=1.9 in /usr/local/lib/python3.6/dist-packages (from cooler) (1.14.6)\n",
            "Requirement already satisfied: pypairix in /usr/local/lib/python3.6/dist-packages (from cooler) (0.3.6)\n",
            "Requirement already satisfied: scipy>=0.16 in /usr/local/lib/python3.6/dist-packages (from cooler) (1.1.0)\n",
            "Requirement already satisfied: h5py>=2.5 in /usr/local/lib/python3.6/dist-packages (from cooler) (2.8.0)\n",
            "Requirement already satisfied: biopython in /usr/local/lib/python3.6/dist-packages (from cooler) (1.72)\n",
            "Requirement already satisfied: six in /usr/local/lib/python3.6/dist-packages (from cooler) (1.11.0)\n",
            "Requirement already satisfied: python-dateutil>=2 in /usr/local/lib/python3.6/dist-packages (from pandas>=0.19->cooler) (2.5.3)\n",
            "Requirement already satisfied: pytz>=2011k in /usr/local/lib/python3.6/dist-packages (from pandas>=0.19->cooler) (2018.7)\n",
            "Requirement already satisfied: setuptools>=0.7 in /usr/local/lib/python3.6/dist-packages (from pyfaidx->cooler) (40.6.2)\n",
            "Requirement already satisfied: dill>=0.2.8.1 in /usr/local/lib/python3.6/dist-packages (from multiprocess->cooler) (0.2.8.2)\n",
            "Requirement already satisfied: toolz>=0.8.0 in /usr/local/lib/python3.6/dist-packages (from cytoolz->cooler) (0.9.0)\n",
            "Requirement already satisfied: bioframe in /usr/local/lib/python3.6/dist-packages (0.0.7)\n",
            "Requirement already satisfied: requests in /usr/local/lib/python3.6/dist-packages (from bioframe) (2.18.4)\n",
            "Requirement already satisfied: pypairix in /usr/local/lib/python3.6/dist-packages (from bioframe) (0.3.6)\n",
            "Requirement already satisfied: pyfaidx in /usr/local/lib/python3.6/dist-packages (from bioframe) (0.5.5.2)\n",
            "Requirement already satisfied: pysam in /usr/local/lib/python3.6/dist-packages (from bioframe) (0.15.1)\n",
            "Requirement already satisfied: six in /usr/local/lib/python3.6/dist-packages (from bioframe) (1.11.0)\n",
            "Requirement already satisfied: pandas>=0.17 in /usr/local/lib/python3.6/dist-packages (from bioframe) (0.22.0)\n",
            "Requirement already satisfied: numpy>=1.9 in /usr/local/lib/python3.6/dist-packages (from bioframe) (1.14.6)\n",
            "Requirement already satisfied: idna<2.7,>=2.5 in /usr/local/lib/python3.6/dist-packages (from requests->bioframe) (2.6)\n",
            "Requirement already satisfied: urllib3<1.23,>=1.21.1 in /usr/local/lib/python3.6/dist-packages (from requests->bioframe) (1.22)\n",
            "Requirement already satisfied: certifi>=2017.4.17 in /usr/local/lib/python3.6/dist-packages (from requests->bioframe) (2018.11.29)\n",
            "Requirement already satisfied: chardet<3.1.0,>=3.0.2 in /usr/local/lib/python3.6/dist-packages (from requests->bioframe) (3.0.4)\n",
            "Requirement already satisfied: setuptools>=0.7 in /usr/local/lib/python3.6/dist-packages (from pyfaidx->bioframe) (40.6.2)\n",
            "Requirement already satisfied: pytz>=2011k in /usr/local/lib/python3.6/dist-packages (from pandas>=0.17->bioframe) (2018.7)\n",
            "Requirement already satisfied: python-dateutil>=2 in /usr/local/lib/python3.6/dist-packages (from pandas>=0.17->bioframe) (2.5.3)\n",
            "Collecting cython\n",
            "\u001b[?25l  Downloading https://files.pythonhosted.org/packages/a2/e0/0592be5b851c8013aa253592606ca265862d27444d908e029fd75d563c9c/Cython-0.29.1-cp36-cp36m-manylinux1_x86_64.whl (2.1MB)\n",
            "\u001b[K    100% |████████████████████████████████| 2.1MB 10.9MB/s \n",
            "\u001b[?25hInstalling collected packages: cython\n",
            "Successfully installed cython-0.29.1\n",
            "Collecting git+https://github.com/mirnylab/cooltools\n",
            "  Cloning https://github.com/mirnylab/cooltools to /tmp/pip-req-build-ew9d1lxa\n",
            "Requirement already satisfied: numpy in /usr/local/lib/python3.6/dist-packages (from cooltools==0.1.0) (1.14.6)\n",
            "Requirement already satisfied: cython in /usr/local/lib/python3.6/dist-packages (from cooltools==0.1.0) (0.29.1)\n",
            "Collecting dask (from cooltools==0.1.0)\n",
            "\u001b[?25l  Downloading https://files.pythonhosted.org/packages/b4/e5/6407c2349622699eab0881e6bc6e978b8da67872a105e5e20b72ff190c65/dask-1.0.0-py2.py3-none-any.whl (685kB)\n",
            "\u001b[K    100% |████████████████████████████████| 686kB 18.6MB/s \n",
            "\u001b[?25hRequirement already satisfied: cytoolz in /usr/local/lib/python3.6/dist-packages (from cooltools==0.1.0) (0.9.0.1)\n",
            "Collecting numba (from cooltools==0.1.0)\n",
            "\u001b[?25l  Downloading https://files.pythonhosted.org/packages/31/55/938f0023a4f37fe24460d46846670aba8170a6b736f1693293e710d4a6d0/numba-0.41.0-cp36-cp36m-manylinux1_x86_64.whl (3.2MB)\n",
            "\u001b[K    100% |████████████████████████████████| 3.2MB 9.0MB/s \n",
            "\u001b[?25hRequirement already satisfied: click>=7 in /usr/local/lib/python3.6/dist-packages (from cooltools==0.1.0) (7.0)\n",
            "Requirement already satisfied: cooler>=0.7 in /usr/local/lib/python3.6/dist-packages (from cooltools==0.1.0) (0.7.11)\n",
            "Requirement already satisfied: bioframe in /usr/local/lib/python3.6/dist-packages (from cooltools==0.1.0) (0.0.7)\n",
            "Requirement already satisfied: scikit-learn in /usr/local/lib/python3.6/dist-packages (from cooltools==0.1.0) (0.19.2)\n",
            "Collecting tables (from cooltools==0.1.0)\n",
            "\u001b[?25l  Downloading https://files.pythonhosted.org/packages/d7/1b/21f4c7f296b718575c17ef25e61c05742a283c45077b4c8d5a190b3e0b59/tables-3.4.4-cp36-cp36m-manylinux1_x86_64.whl (3.8MB)\n",
            "\u001b[K    100% |████████████████████████████████| 3.8MB 10.6MB/s \n",
            "\u001b[?25hRequirement already satisfied: toolz>=0.8.0 in /usr/local/lib/python3.6/dist-packages (from cytoolz->cooltools==0.1.0) (0.9.0)\n",
            "Collecting llvmlite>=0.26.0dev0 (from numba->cooltools==0.1.0)\n",
            "\u001b[?25l  Downloading https://files.pythonhosted.org/packages/a2/60/d22966c97a47687ac1cc57c2e756380897c264f1ce40780105d7dbcd9564/llvmlite-0.26.0-cp36-cp36m-manylinux1_x86_64.whl (16.1MB)\n",
            "\u001b[K    100% |████████████████████████████████| 16.1MB 2.1MB/s \n",
            "\u001b[?25hRequirement already satisfied: pandas>=0.19 in /usr/local/lib/python3.6/dist-packages (from cooler>=0.7->cooltools==0.1.0) (0.22.0)\n",
            "Requirement already satisfied: biopython in /usr/local/lib/python3.6/dist-packages (from cooler>=0.7->cooltools==0.1.0) (1.72)\n",
            "Requirement already satisfied: six in /usr/local/lib/python3.6/dist-packages (from cooler>=0.7->cooltools==0.1.0) (1.11.0)\n",
            "Requirement already satisfied: pypairix in /usr/local/lib/python3.6/dist-packages (from cooler>=0.7->cooltools==0.1.0) (0.3.6)\n",
            "Requirement already satisfied: multiprocess in /usr/local/lib/python3.6/dist-packages (from cooler>=0.7->cooltools==0.1.0) (0.70.6.1)\n",
            "Requirement already satisfied: scipy>=0.16 in /usr/local/lib/python3.6/dist-packages (from cooler>=0.7->cooltools==0.1.0) (1.1.0)\n",
            "Requirement already satisfied: pyfaidx in /usr/local/lib/python3.6/dist-packages (from cooler>=0.7->cooltools==0.1.0) (0.5.5.2)\n",
            "Requirement already satisfied: pysam>0.8 in /usr/local/lib/python3.6/dist-packages (from cooler>=0.7->cooltools==0.1.0) (0.15.1)\n",
            "Requirement already satisfied: h5py>=2.5 in /usr/local/lib/python3.6/dist-packages (from cooler>=0.7->cooltools==0.1.0) (2.8.0)\n",
            "Requirement already satisfied: requests in /usr/local/lib/python3.6/dist-packages (from bioframe->cooltools==0.1.0) (2.18.4)\n",
            "Collecting numexpr>=2.5.2 (from tables->cooltools==0.1.0)\n",
            "\u001b[?25l  Downloading https://files.pythonhosted.org/packages/db/ea/efd9e16283637eb5b6c0042b6cc3521f1b9a5b47767ac463c88bbd37670c/numexpr-2.6.8-cp36-cp36m-manylinux1_x86_64.whl (162kB)\n",
            "\u001b[K    100% |████████████████████████████████| 163kB 24.2MB/s \n",
            "\u001b[?25hRequirement already satisfied: pytz>=2011k in /usr/local/lib/python3.6/dist-packages (from pandas>=0.19->cooler>=0.7->cooltools==0.1.0) (2018.7)\n",
            "Requirement already satisfied: python-dateutil>=2 in /usr/local/lib/python3.6/dist-packages (from pandas>=0.19->cooler>=0.7->cooltools==0.1.0) (2.5.3)\n",
            "Requirement already satisfied: dill>=0.2.8.1 in /usr/local/lib/python3.6/dist-packages (from multiprocess->cooler>=0.7->cooltools==0.1.0) (0.2.8.2)\n",
            "Requirement already satisfied: setuptools>=0.7 in /usr/local/lib/python3.6/dist-packages (from pyfaidx->cooler>=0.7->cooltools==0.1.0) (40.6.2)\n",
            "Requirement already satisfied: certifi>=2017.4.17 in /usr/local/lib/python3.6/dist-packages (from requests->bioframe->cooltools==0.1.0) (2018.11.29)\n",
            "Requirement already satisfied: urllib3<1.23,>=1.21.1 in /usr/local/lib/python3.6/dist-packages (from requests->bioframe->cooltools==0.1.0) (1.22)\n",
            "Requirement already satisfied: idna<2.7,>=2.5 in /usr/local/lib/python3.6/dist-packages (from requests->bioframe->cooltools==0.1.0) (2.6)\n",
            "Requirement already satisfied: chardet<3.1.0,>=3.0.2 in /usr/local/lib/python3.6/dist-packages (from requests->bioframe->cooltools==0.1.0) (3.0.4)\n",
            "Building wheels for collected packages: cooltools\n",
            "  Running setup.py bdist_wheel for cooltools ... \u001b[?25l-\b \b\\\b \b|\b \b/\b \b-\b \b\\\b \b|\b \bdone\n",
            "\u001b[?25h  Stored in directory: /tmp/pip-ephem-wheel-cache-89z1dn1v/wheels/bf/59/ac/f23113200c2678fc1dd5901a905c180fc6392ae3af501292d0\n",
            "Successfully built cooltools\n",
            "Installing collected packages: dask, llvmlite, numba, numexpr, tables, cooltools\n",
            "Successfully installed cooltools-0.1.0 dask-1.0.0 llvmlite-0.26.0 numba-0.41.0 numexpr-2.6.8 tables-3.4.4\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "metadata": {
        "id": "fH3ryJV7vpFE",
        "colab_type": "code",
        "outputId": "a9664612-edd0-4b0d-e6c7-3439964c44e2",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 204
        }
      },
      "cell_type": "code",
      "source": [
        "!wget http://hgdownload.soe.ucsc.edu/admin/exe/linux.x86_64.v369/bedGraphToBigWig\n",
        "!chmod 755 bedGraphToBigWig"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "--2018-12-11 19:45:38--  http://hgdownload.soe.ucsc.edu/admin/exe/linux.x86_64.v369/bedGraphToBigWig\n",
            "Resolving hgdownload.soe.ucsc.edu (hgdownload.soe.ucsc.edu)... 128.114.119.163\n",
            "Connecting to hgdownload.soe.ucsc.edu (hgdownload.soe.ucsc.edu)|128.114.119.163|:80... connected.\n",
            "HTTP request sent, awaiting response... 200 OK\n",
            "Length: 4380534 (4.2M) [text/plain]\n",
            "Saving to: ‘bedGraphToBigWig.2’\n",
            "\n",
            "bedGraphToBigWig.2  100%[===================>]   4.18M  3.22MB/s    in 1.3s    \n",
            "\n",
            "2018-12-11 19:45:39 (3.22 MB/s) - ‘bedGraphToBigWig.2’ saved [4380534/4380534]\n",
            "\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "metadata": {
        "id": "H92Qjhb17WqU",
        "colab_type": "code",
        "colab": {}
      },
      "cell_type": "code",
      "source": [
        "cp bedGraphToBigWig /usr/bin"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "metadata": {
        "id": "eUJ0T0ow7DGi",
        "colab_type": "code",
        "outputId": "69e7fa3b-2c44-4c59-e629-44697e830c18",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 323
        }
      },
      "cell_type": "code",
      "source": [
        "#Note: !bedGraphToBigWig needs to be executable (i.e. not !./bedGraphToBigWig)\n",
        "!bedGraphToBigWig "
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "bedGraphToBigWig v 4 - Convert a bedGraph file to bigWig format.\n",
            "usage:\n",
            "   bedGraphToBigWig in.bedGraph chrom.sizes out.bw\n",
            "where in.bedGraph is a four column file in the format:\n",
            "      <chrom> <start> <end> <value>\n",
            "and chrom.sizes is a two-column file/URL: <chromosome name> <size in bases>\n",
            "and out.bw is the output indexed big wig file.\n",
            "If the assembly <db> is hosted by UCSC, chrom.sizes can be a URL like\n",
            "  http://hgdownload.cse.ucsc.edu/goldenPath/<db>/bigZips/<db>.chrom.sizes\n",
            "or you may use the script fetchChromSizes to download the chrom.sizes file.\n",
            "If not hosted by UCSC, a chrom.sizes file can be generated by running\n",
            "twoBitInfo on the assembly .2bit file.\n",
            "The input bedGraph file must be sorted, use the unix sort command:\n",
            "  sort -k1,1 -k2,2n unsorted.bedGraph > sorted.bedGraph\n",
            "options:\n",
            "   -blockSize=N - Number of items to bundle in r-tree.  Default 256\n",
            "   -itemsPerSlot=N - Number of data points bundled at lowest level. Default 1024\n",
            "   -unc - If set, do not use compression.\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "metadata": {
        "id": "_qdppq0W_yEc",
        "colab_type": "text"
      },
      "cell_type": "markdown",
      "source": [
        "###Upload .cool  Files:"
      ]
    },
    {
      "metadata": {
        "id": "hZa4Ugsx503B",
        "colab_type": "code",
        "outputId": "c5990129-2298-4ad6-a17c-f6bad1338bd6",
        "colab": {
          "resources": {
            "http://localhost:8080/nbextensions/google.colab/files.js": {
              "data": "Ly8gQ29weXJpZ2h0IDIwMTcgR29vZ2xlIExMQwovLwovLyBMaWNlbnNlZCB1bmRlciB0aGUgQXBhY2hlIExpY2Vuc2UsIFZlcnNpb24gMi4wICh0aGUgIkxpY2Vuc2UiKTsKLy8geW91IG1heSBub3QgdXNlIHRoaXMgZmlsZSBleGNlcHQgaW4gY29tcGxpYW5jZSB3aXRoIHRoZSBMaWNlbnNlLgovLyBZb3UgbWF5IG9idGFpbiBhIGNvcHkgb2YgdGhlIExpY2Vuc2UgYXQKLy8KLy8gICAgICBodHRwOi8vd3d3LmFwYWNoZS5vcmcvbGljZW5zZXMvTElDRU5TRS0yLjAKLy8KLy8gVW5sZXNzIHJlcXVpcmVkIGJ5IGFwcGxpY2FibGUgbGF3IG9yIGFncmVlZCB0byBpbiB3cml0aW5nLCBzb2Z0d2FyZQovLyBkaXN0cmlidXRlZCB1bmRlciB0aGUgTGljZW5zZSBpcyBkaXN0cmlidXRlZCBvbiBhbiAiQVMgSVMiIEJBU0lTLAovLyBXSVRIT1VUIFdBUlJBTlRJRVMgT1IgQ09ORElUSU9OUyBPRiBBTlkgS0lORCwgZWl0aGVyIGV4cHJlc3Mgb3IgaW1wbGllZC4KLy8gU2VlIHRoZSBMaWNlbnNlIGZvciB0aGUgc3BlY2lmaWMgbGFuZ3VhZ2UgZ292ZXJuaW5nIHBlcm1pc3Npb25zIGFuZAovLyBsaW1pdGF0aW9ucyB1bmRlciB0aGUgTGljZW5zZS4KCi8qKgogKiBAZmlsZW92ZXJ2aWV3IEhlbHBlcnMgZm9yIGdvb2dsZS5jb2xhYiBQeXRob24gbW9kdWxlLgogKi8KKGZ1bmN0aW9uKHNjb3BlKSB7CmZ1bmN0aW9uIHNwYW4odGV4dCwgc3R5bGVBdHRyaWJ1dGVzID0ge30pIHsKICBjb25zdCBlbGVtZW50ID0gZG9jdW1lbnQuY3JlYXRlRWxlbWVudCgnc3BhbicpOwogIGVsZW1lbnQudGV4dENvbnRlbnQgPSB0ZXh0OwogIGZvciAoY29uc3Qga2V5IG9mIE9iamVjdC5rZXlzKHN0eWxlQXR0cmlidXRlcykpIHsKICAgIGVsZW1lbnQuc3R5bGVba2V5XSA9IHN0eWxlQXR0cmlidXRlc1trZXldOwogIH0KICByZXR1cm4gZWxlbWVudDsKfQoKLy8gTWF4IG51bWJlciBvZiBieXRlcyB3aGljaCB3aWxsIGJlIHVwbG9hZGVkIGF0IGEgdGltZS4KY29uc3QgTUFYX1BBWUxPQURfU0laRSA9IDEwMCAqIDEwMjQ7Ci8vIE1heCBhbW91bnQgb2YgdGltZSB0byBibG9jayB3YWl0aW5nIGZvciB0aGUgdXNlci4KY29uc3QgRklMRV9DSEFOR0VfVElNRU9VVF9NUyA9IDMwICogMTAwMDsKCmZ1bmN0aW9uIF91cGxvYWRGaWxlcyhpbnB1dElkLCBvdXRwdXRJZCkgewogIGNvbnN0IHN0ZXBzID0gdXBsb2FkRmlsZXNTdGVwKGlucHV0SWQsIG91dHB1dElkKTsKICBjb25zdCBvdXRwdXRFbGVtZW50ID0gZG9jdW1lbnQuZ2V0RWxlbWVudEJ5SWQob3V0cHV0SWQpOwogIC8vIENhY2hlIHN0ZXBzIG9uIHRoZSBvdXRwdXRFbGVtZW50IHRvIG1ha2UgaXQgYXZhaWxhYmxlIGZvciB0aGUgbmV4dCBjYWxsCiAgLy8gdG8gdXBsb2FkRmlsZXNDb250aW51ZSBmcm9tIFB5dGhvbi4KICBvdXRwdXRFbGVtZW50LnN0ZXBzID0gc3RlcHM7CgogIHJldHVybiBfdXBsb2FkRmlsZXNDb250aW51ZShvdXRwdXRJZCk7Cn0KCi8vIFRoaXMgaXMgcm91Z2hseSBhbiBhc3luYyBnZW5lcmF0b3IgKG5vdCBzdXBwb3J0ZWQgaW4gdGhlIGJyb3dzZXIgeWV0KSwKLy8gd2hlcmUgdGhlcmUgYXJlIG11bHRpcGxlIGFzeW5jaHJvbm91cyBzdGVwcyBhbmQgdGhlIFB5dGhvbiBzaWRlIGlzIGdvaW5nCi8vIHRvIHBvbGwgZm9yIGNvbXBsZXRpb24gb2YgZWFjaCBzdGVwLgovLyBUaGlzIHVzZXMgYSBQcm9taXNlIHRvIGJsb2NrIHRoZSBweXRob24gc2lkZSBvbiBjb21wbGV0aW9uIG9mIGVhY2ggc3RlcCwKLy8gdGhlbiBwYXNzZXMgdGhlIHJlc3VsdCBvZiB0aGUgcHJldmlvdXMgc3RlcCBhcyB0aGUgaW5wdXQgdG8gdGhlIG5leHQgc3RlcC4KZnVuY3Rpb24gX3VwbG9hZEZpbGVzQ29udGludWUob3V0cHV0SWQpIHsKICBjb25zdCBvdXRwdXRFbGVtZW50ID0gZG9jdW1lbnQuZ2V0RWxlbWVudEJ5SWQob3V0cHV0SWQpOwogIGNvbnN0IHN0ZXBzID0gb3V0cHV0RWxlbWVudC5zdGVwczsKCiAgY29uc3QgbmV4dCA9IHN0ZXBzLm5leHQob3V0cHV0RWxlbWVudC5sYXN0UHJvbWlzZVZhbHVlKTsKICByZXR1cm4gUHJvbWlzZS5yZXNvbHZlKG5leHQudmFsdWUucHJvbWlzZSkudGhlbigodmFsdWUpID0+IHsKICAgIC8vIENhY2hlIHRoZSBsYXN0IHByb21pc2UgdmFsdWUgdG8gbWFrZSBpdCBhdmFpbGFibGUgdG8gdGhlIG5leHQKICAgIC8vIHN0ZXAgb2YgdGhlIGdlbmVyYXRvci4KICAgIG91dHB1dEVsZW1lbnQubGFzdFByb21pc2VWYWx1ZSA9IHZhbHVlOwogICAgcmV0dXJuIG5leHQudmFsdWUucmVzcG9uc2U7CiAgfSk7Cn0KCi8qKgogKiBHZW5lcmF0b3IgZnVuY3Rpb24gd2hpY2ggaXMgY2FsbGVkIGJldHdlZW4gZWFjaCBhc3luYyBzdGVwIG9mIHRoZSB1cGxvYWQKICogcHJvY2Vzcy4KICogQHBhcmFtIHtzdHJpbmd9IGlucHV0SWQgRWxlbWVudCBJRCBvZiB0aGUgaW5wdXQgZmlsZSBwaWNrZXIgZWxlbWVudC4KICogQHBhcmFtIHtzdHJpbmd9IG91dHB1dElkIEVsZW1lbnQgSUQgb2YgdGhlIG91dHB1dCBkaXNwbGF5LgogKiBAcmV0dXJuIHshSXRlcmFibGU8IU9iamVjdD59IEl0ZXJhYmxlIG9mIG5leHQgc3RlcHMuCiAqLwpmdW5jdGlvbiogdXBsb2FkRmlsZXNTdGVwKGlucHV0SWQsIG91dHB1dElkKSB7CiAgY29uc3QgaW5wdXRFbGVtZW50ID0gZG9jdW1lbnQuZ2V0RWxlbWVudEJ5SWQoaW5wdXRJZCk7CiAgaW5wdXRFbGVtZW50LmRpc2FibGVkID0gZmFsc2U7CgogIGNvbnN0IG91dHB1dEVsZW1lbnQgPSBkb2N1bWVudC5nZXRFbGVtZW50QnlJZChvdXRwdXRJZCk7CiAgb3V0cHV0RWxlbWVudC5pbm5lckhUTUwgPSAnJzsKCiAgY29uc3QgcGlja2VkUHJvbWlzZSA9IG5ldyBQcm9taXNlKChyZXNvbHZlKSA9PiB7CiAgICBpbnB1dEVsZW1lbnQuYWRkRXZlbnRMaXN0ZW5lcignY2hhbmdlJywgKGUpID0+IHsKICAgICAgcmVzb2x2ZShlLnRhcmdldC5maWxlcyk7CiAgICB9KTsKICB9KTsKCiAgY29uc3QgY2FuY2VsID0gZG9jdW1lbnQuY3JlYXRlRWxlbWVudCgnYnV0dG9uJyk7CiAgaW5wdXRFbGVtZW50LnBhcmVudEVsZW1lbnQuYXBwZW5kQ2hpbGQoY2FuY2VsKTsKICBjYW5jZWwudGV4dENvbnRlbnQgPSAnQ2FuY2VsIHVwbG9hZCc7CiAgY29uc3QgY2FuY2VsUHJvbWlzZSA9IG5ldyBQcm9taXNlKChyZXNvbHZlKSA9PiB7CiAgICBjYW5jZWwub25jbGljayA9ICgpID0+IHsKICAgICAgcmVzb2x2ZShudWxsKTsKICAgIH07CiAgfSk7CgogIC8vIENhbmNlbCB1cGxvYWQgaWYgdXNlciBoYXNuJ3QgcGlja2VkIGFueXRoaW5nIGluIHRpbWVvdXQuCiAgY29uc3QgdGltZW91dFByb21pc2UgPSBuZXcgUHJvbWlzZSgocmVzb2x2ZSkgPT4gewogICAgc2V0VGltZW91dCgoKSA9PiB7CiAgICAgIHJlc29sdmUobnVsbCk7CiAgICB9LCBGSUxFX0NIQU5HRV9USU1FT1VUX01TKTsKICB9KTsKCiAgLy8gV2FpdCBmb3IgdGhlIHVzZXIgdG8gcGljayB0aGUgZmlsZXMuCiAgY29uc3QgZmlsZXMgPSB5aWVsZCB7CiAgICBwcm9taXNlOiBQcm9taXNlLnJhY2UoW3BpY2tlZFByb21pc2UsIHRpbWVvdXRQcm9taXNlLCBjYW5jZWxQcm9taXNlXSksCiAgICByZXNwb25zZTogewogICAgICBhY3Rpb246ICdzdGFydGluZycsCiAgICB9CiAgfTsKCiAgaWYgKCFmaWxlcykgewogICAgcmV0dXJuIHsKICAgICAgcmVzcG9uc2U6IHsKICAgICAgICBhY3Rpb246ICdjb21wbGV0ZScsCiAgICAgIH0KICAgIH07CiAgfQoKICBjYW5jZWwucmVtb3ZlKCk7CgogIC8vIERpc2FibGUgdGhlIGlucHV0IGVsZW1lbnQgc2luY2UgZnVydGhlciBwaWNrcyBhcmUgbm90IGFsbG93ZWQuCiAgaW5wdXRFbGVtZW50LmRpc2FibGVkID0gdHJ1ZTsKCiAgZm9yIChjb25zdCBmaWxlIG9mIGZpbGVzKSB7CiAgICBjb25zdCBsaSA9IGRvY3VtZW50LmNyZWF0ZUVsZW1lbnQoJ2xpJyk7CiAgICBsaS5hcHBlbmQoc3BhbihmaWxlLm5hbWUsIHtmb250V2VpZ2h0OiAnYm9sZCd9KSk7CiAgICBsaS5hcHBlbmQoc3BhbigKICAgICAgICBgKCR7ZmlsZS50eXBlIHx8ICduL2EnfSkgLSAke2ZpbGUuc2l6ZX0gYnl0ZXMsIGAgKwogICAgICAgIGBsYXN0IG1vZGlmaWVkOiAkewogICAgICAgICAgICBmaWxlLmxhc3RNb2RpZmllZERhdGUgPyBmaWxlLmxhc3RNb2RpZmllZERhdGUudG9Mb2NhbGVEYXRlU3RyaW5nKCkgOgogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAnbi9hJ30gLSBgKSk7CiAgICBjb25zdCBwZXJjZW50ID0gc3BhbignMCUgZG9uZScpOwogICAgbGkuYXBwZW5kQ2hpbGQocGVyY2VudCk7CgogICAgb3V0cHV0RWxlbWVudC5hcHBlbmRDaGlsZChsaSk7CgogICAgY29uc3QgZmlsZURhdGFQcm9taXNlID0gbmV3IFByb21pc2UoKHJlc29sdmUpID0+IHsKICAgICAgY29uc3QgcmVhZGVyID0gbmV3IEZpbGVSZWFkZXIoKTsKICAgICAgcmVhZGVyLm9ubG9hZCA9IChlKSA9PiB7CiAgICAgICAgcmVzb2x2ZShlLnRhcmdldC5yZXN1bHQpOwogICAgICB9OwogICAgICByZWFkZXIucmVhZEFzQXJyYXlCdWZmZXIoZmlsZSk7CiAgICB9KTsKICAgIC8vIFdhaXQgZm9yIHRoZSBkYXRhIHRvIGJlIHJlYWR5LgogICAgbGV0IGZpbGVEYXRhID0geWllbGQgewogICAgICBwcm9taXNlOiBmaWxlRGF0YVByb21pc2UsCiAgICAgIHJlc3BvbnNlOiB7CiAgICAgICAgYWN0aW9uOiAnY29udGludWUnLAogICAgICB9CiAgICB9OwoKICAgIC8vIFVzZSBhIGNodW5rZWQgc2VuZGluZyB0byBhdm9pZCBtZXNzYWdlIHNpemUgbGltaXRzLiBTZWUgYi82MjExNTY2MC4KICAgIGxldCBwb3NpdGlvbiA9IDA7CiAgICB3aGlsZSAocG9zaXRpb24gPCBmaWxlRGF0YS5ieXRlTGVuZ3RoKSB7CiAgICAgIGNvbnN0IGxlbmd0aCA9IE1hdGgubWluKGZpbGVEYXRhLmJ5dGVMZW5ndGggLSBwb3NpdGlvbiwgTUFYX1BBWUxPQURfU0laRSk7CiAgICAgIGNvbnN0IGNodW5rID0gbmV3IFVpbnQ4QXJyYXkoZmlsZURhdGEsIHBvc2l0aW9uLCBsZW5ndGgpOwogICAgICBwb3NpdGlvbiArPSBsZW5ndGg7CgogICAgICBjb25zdCBiYXNlNjQgPSBidG9hKFN0cmluZy5mcm9tQ2hhckNvZGUuYXBwbHkobnVsbCwgY2h1bmspKTsKICAgICAgeWllbGQgewogICAgICAgIHJlc3BvbnNlOiB7CiAgICAgICAgICBhY3Rpb246ICdhcHBlbmQnLAogICAgICAgICAgZmlsZTogZmlsZS5uYW1lLAogICAgICAgICAgZGF0YTogYmFzZTY0LAogICAgICAgIH0sCiAgICAgIH07CiAgICAgIHBlcmNlbnQudGV4dENvbnRlbnQgPQogICAgICAgICAgYCR7TWF0aC5yb3VuZCgocG9zaXRpb24gLyBmaWxlRGF0YS5ieXRlTGVuZ3RoKSAqIDEwMCl9JSBkb25lYDsKICAgIH0KICB9CgogIC8vIEFsbCBkb25lLgogIHlpZWxkIHsKICAgIHJlc3BvbnNlOiB7CiAgICAgIGFjdGlvbjogJ2NvbXBsZXRlJywKICAgIH0KICB9Owp9CgpzY29wZS5nb29nbGUgPSBzY29wZS5nb29nbGUgfHwge307CnNjb3BlLmdvb2dsZS5jb2xhYiA9IHNjb3BlLmdvb2dsZS5jb2xhYiB8fCB7fTsKc2NvcGUuZ29vZ2xlLmNvbGFiLl9maWxlcyA9IHsKICBfdXBsb2FkRmlsZXMsCiAgX3VwbG9hZEZpbGVzQ29udGludWUsCn07Cn0pKHNlbGYpOwo=",
              "ok": true,
              "headers": [
                [
                  "content-type",
                  "application/javascript"
                ]
              ],
              "status": 200,
              "status_text": ""
            }
          },
          "base_uri": "https://localhost:8080/",
          "height": 35
        }
      },
      "cell_type": "code",
      "source": [
        "from google.colab import files\n",
        "\n",
        "uploaded = files.upload()\n",
        "\n",
        "for fn in uploaded.keys():\n",
        "  print('User uploaded file \"{name}\" with length {length} bytes'.format(name=fn, length=len(uploaded[fn])))"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "text/html": [
              "\n",
              "     <input type=\"file\" id=\"files-b15477f5-1fa9-4583-a1b7-b1f128b477ce\" name=\"files[]\" multiple disabled />\n",
              "     <output id=\"result-b15477f5-1fa9-4583-a1b7-b1f128b477ce\">\n",
              "      Upload widget is only available when the cell has been executed in the\n",
              "      current browser session. Please rerun this cell to enable.\n",
              "      </output>\n",
              "      <script src=\"/nbextensions/google.colab/files.js\"></script> "
            ],
            "text/plain": [
              "<IPython.core.display.HTML object>"
            ]
          },
          "metadata": {
            "tags": []
          }
        }
      ]
    },
    {
      "metadata": {
        "id": "ZoLg_--RvrMc",
        "colab_type": "text"
      },
      "cell_type": "markdown",
      "source": [
        "### chromarm eigenvector"
      ]
    },
    {
      "metadata": {
        "id": "niH0IK_vu2h1",
        "colab_type": "text"
      },
      "cell_type": "markdown",
      "source": [
        "Steps:\n",
        "\n",
        "* Download `twoBitToFa` and the hg19 sequence\n",
        "* Convert to fasta\n",
        "* Download hg19 gapfile\n",
        "* Generate chromarm table to get the `regions`\n",
        "* Get the bin table from the cooler\n",
        "* Load the fasta file and calculate GC content on your bin table\n",
        "* Pass the bin table with GC and the list of `regions` to `cooltools.eigdecomp.cooler_cis_eig`"
      ]
    },
    {
      "metadata": {
        "id": "lr33-1zgtnky",
        "colab_type": "code",
        "outputId": "37d0a3d5-db11-473e-8676-2b2c9baaf174",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 204
        }
      },
      "cell_type": "code",
      "source": [
        "!wget http://hgdownload.soe.ucsc.edu/admin/exe/linux.x86_64.v369/twoBitToFa\n",
        "!chmod 755 twoBitToFa"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "--2018-12-11 19:45:42--  http://hgdownload.soe.ucsc.edu/admin/exe/linux.x86_64.v369/twoBitToFa\n",
            "Resolving hgdownload.soe.ucsc.edu (hgdownload.soe.ucsc.edu)... 128.114.119.163\n",
            "Connecting to hgdownload.soe.ucsc.edu (hgdownload.soe.ucsc.edu)|128.114.119.163|:80... connected.\n",
            "HTTP request sent, awaiting response... 200 OK\n",
            "Length: 4716626 (4.5M) [text/plain]\n",
            "Saving to: ‘twoBitToFa’\n",
            "\n",
            "twoBitToFa          100%[===================>]   4.50M  3.58MB/s    in 1.3s    \n",
            "\n",
            "2018-12-11 19:45:43 (3.58 MB/s) - ‘twoBitToFa’ saved [4716626/4716626]\n",
            "\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "metadata": {
        "id": "ANaMINeFsxWC",
        "colab_type": "code",
        "outputId": "8ef80b05-f353-4f19-b4b5-003738ca1182",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 204
        }
      },
      "cell_type": "code",
      "source": [
        "!wget http://hgdownload.cse.ucsc.edu/goldenPath/hg19/bigZips/hg19.2bit"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "--2018-12-11 19:45:46--  http://hgdownload.cse.ucsc.edu/goldenPath/hg19/bigZips/hg19.2bit\n",
            "Resolving hgdownload.cse.ucsc.edu (hgdownload.cse.ucsc.edu)... 128.114.119.163\n",
            "Connecting to hgdownload.cse.ucsc.edu (hgdownload.cse.ucsc.edu)|128.114.119.163|:80... connected.\n",
            "HTTP request sent, awaiting response... 200 OK\n",
            "Length: 816241703 (778M) [text/plain]\n",
            "Saving to: ‘hg19.2bit’\n",
            "\n",
            "hg19.2bit           100%[===================>] 778.43M  30.9MB/s    in 16s     \n",
            "\n",
            "2018-12-11 19:46:03 (47.7 MB/s) - ‘hg19.2bit’ saved [816241703/816241703]\n",
            "\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "metadata": {
        "id": "VIuAM2Hjtq9h",
        "colab_type": "code",
        "colab": {}
      },
      "cell_type": "code",
      "source": [
        "!./twoBitToFa hg19.2bit hg19.fa"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "metadata": {
        "id": "UvNcYWeZq65X",
        "colab_type": "code",
        "colab": {}
      },
      "cell_type": "code",
      "source": [
        "import bioframe\n",
        "import cooltools.eigdecomp"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "metadata": {
        "id": "OJqnmxLfqdnd",
        "colab_type": "code",
        "outputId": "d0bc6e54-a784-4120-b40b-0a9c90285908",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 204
        }
      },
      "cell_type": "code",
      "source": [
        "hg19 = bioframe.fetch_chromsizes('hg19')\n",
        "gaps = bioframe.read_gapfile('http://hgdownload.cse.ucsc.edu/goldenPath/hg19/database/gap.txt.gz')\n",
        "gaps = gaps[gaps['type'] == 'centromere']\n",
        "gaps['mid'] = (gaps['start'] + gaps['end']) // 2\n",
        "gaps.head()"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/html": [
              "<div>\n",
              "<style scoped>\n",
              "    .dataframe tbody tr th:only-of-type {\n",
              "        vertical-align: middle;\n",
              "    }\n",
              "\n",
              "    .dataframe tbody tr th {\n",
              "        vertical-align: top;\n",
              "    }\n",
              "\n",
              "    .dataframe thead th {\n",
              "        text-align: right;\n",
              "    }\n",
              "</style>\n",
              "<table border=\"1\" class=\"dataframe\">\n",
              "  <thead>\n",
              "    <tr style=\"text-align: right;\">\n",
              "      <th></th>\n",
              "      <th>chrom</th>\n",
              "      <th>start</th>\n",
              "      <th>end</th>\n",
              "      <th>length</th>\n",
              "      <th>type</th>\n",
              "      <th>bridge</th>\n",
              "      <th>mid</th>\n",
              "    </tr>\n",
              "  </thead>\n",
              "  <tbody>\n",
              "    <tr>\n",
              "      <th>15</th>\n",
              "      <td>chr1</td>\n",
              "      <td>121535434</td>\n",
              "      <td>124535434</td>\n",
              "      <td>3000000</td>\n",
              "      <td>centromere</td>\n",
              "      <td>no</td>\n",
              "      <td>123035434</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>46</th>\n",
              "      <td>chr10</td>\n",
              "      <td>39254935</td>\n",
              "      <td>42254935</td>\n",
              "      <td>3000000</td>\n",
              "      <td>centromere</td>\n",
              "      <td>no</td>\n",
              "      <td>40754935</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>68</th>\n",
              "      <td>chr11</td>\n",
              "      <td>51644205</td>\n",
              "      <td>54644205</td>\n",
              "      <td>3000000</td>\n",
              "      <td>centromere</td>\n",
              "      <td>no</td>\n",
              "      <td>53144205</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>80</th>\n",
              "      <td>chr12</td>\n",
              "      <td>34856694</td>\n",
              "      <td>37856694</td>\n",
              "      <td>3000000</td>\n",
              "      <td>centromere</td>\n",
              "      <td>no</td>\n",
              "      <td>36356694</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>87</th>\n",
              "      <td>chr13</td>\n",
              "      <td>16000000</td>\n",
              "      <td>19000000</td>\n",
              "      <td>3000000</td>\n",
              "      <td>centromere</td>\n",
              "      <td>no</td>\n",
              "      <td>17500000</td>\n",
              "    </tr>\n",
              "  </tbody>\n",
              "</table>\n",
              "</div>"
            ],
            "text/plain": [
              "    chrom      start        end   length        type bridge        mid\n",
              "15   chr1  121535434  124535434  3000000  centromere     no  123035434\n",
              "46  chr10   39254935   42254935  3000000  centromere     no   40754935\n",
              "68  chr11   51644205   54644205  3000000  centromere     no   53144205\n",
              "80  chr12   34856694   37856694  3000000  centromere     no   36356694\n",
              "87  chr13   16000000   19000000  3000000  centromere     no   17500000"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 45
        }
      ]
    },
    {
      "metadata": {
        "id": "yfG_VrcFqXle",
        "colab_type": "code",
        "outputId": "b0b13ad4-7837-4141-c3a1-52e25b7533ad",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 204
        }
      },
      "cell_type": "code",
      "source": [
        "chromarms = bioframe.tools.split_chromosomes(hg19, gaps.set_index('chrom')['mid'], suffixes=('p', 'q'))\n",
        "chromarms.head()"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/html": [
              "<div>\n",
              "<style scoped>\n",
              "    .dataframe tbody tr th:only-of-type {\n",
              "        vertical-align: middle;\n",
              "    }\n",
              "\n",
              "    .dataframe tbody tr th {\n",
              "        vertical-align: top;\n",
              "    }\n",
              "\n",
              "    .dataframe thead th {\n",
              "        text-align: right;\n",
              "    }\n",
              "</style>\n",
              "<table border=\"1\" class=\"dataframe\">\n",
              "  <thead>\n",
              "    <tr style=\"text-align: right;\">\n",
              "      <th></th>\n",
              "      <th>chrom</th>\n",
              "      <th>start</th>\n",
              "      <th>end</th>\n",
              "      <th>name</th>\n",
              "    </tr>\n",
              "  </thead>\n",
              "  <tbody>\n",
              "    <tr>\n",
              "      <th>0</th>\n",
              "      <td>chr1</td>\n",
              "      <td>0</td>\n",
              "      <td>123035434</td>\n",
              "      <td>chr1p</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>1</th>\n",
              "      <td>chr1</td>\n",
              "      <td>123035434</td>\n",
              "      <td>249250621</td>\n",
              "      <td>chr1q</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>2</th>\n",
              "      <td>chr2</td>\n",
              "      <td>0</td>\n",
              "      <td>93826171</td>\n",
              "      <td>chr2p</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>3</th>\n",
              "      <td>chr2</td>\n",
              "      <td>93826171</td>\n",
              "      <td>243199373</td>\n",
              "      <td>chr2q</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>4</th>\n",
              "      <td>chr3</td>\n",
              "      <td>0</td>\n",
              "      <td>92004854</td>\n",
              "      <td>chr3p</td>\n",
              "    </tr>\n",
              "  </tbody>\n",
              "</table>\n",
              "</div>"
            ],
            "text/plain": [
              "  chrom      start        end   name\n",
              "0  chr1          0  123035434  chr1p\n",
              "1  chr1  123035434  249250621  chr1q\n",
              "2  chr2          0   93826171  chr2p\n",
              "3  chr2   93826171  243199373  chr2q\n",
              "4  chr3          0   92004854  chr3p"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 46
        }
      ]
    },
    {
      "metadata": {
        "id": "dWVYy13G6PbR",
        "colab_type": "code",
        "outputId": "e513294c-e6fc-4c85-e523-28ccea89555a",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 34
        }
      },
      "cell_type": "code",
      "source": [
        "import cooler\n",
        "clr = cooler.Cooler('DKO.20000.cool')\n",
        "clr"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "<Cooler \"DKO.20000.cool::/\">"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 47
        }
      ]
    },
    {
      "metadata": {
        "id": "lwMCB6aey_nk",
        "colab_type": "code",
        "colab": {}
      },
      "cell_type": "code",
      "source": [
        "fasta_records = bioframe.load_fasta('hg19.fa')"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "metadata": {
        "id": "YBQUmT47r0K_",
        "colab_type": "code",
        "outputId": "3b43ebf1-8caa-4e64-ca35-85ca4d7ad4b8",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 204
        }
      },
      "cell_type": "code",
      "source": [
        "bins = clr.bins()[:]\n",
        "fasta_records = bioframe.load_fasta('hg19.fa')\n",
        "bins['GC'] = bioframe.tools.frac_gc(bins, fasta_records)\n",
        "bins.head()"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/html": [
              "<div>\n",
              "<style scoped>\n",
              "    .dataframe tbody tr th:only-of-type {\n",
              "        vertical-align: middle;\n",
              "    }\n",
              "\n",
              "    .dataframe tbody tr th {\n",
              "        vertical-align: top;\n",
              "    }\n",
              "\n",
              "    .dataframe thead th {\n",
              "        text-align: right;\n",
              "    }\n",
              "</style>\n",
              "<table border=\"1\" class=\"dataframe\">\n",
              "  <thead>\n",
              "    <tr style=\"text-align: right;\">\n",
              "      <th></th>\n",
              "      <th>chrom</th>\n",
              "      <th>start</th>\n",
              "      <th>end</th>\n",
              "      <th>weight</th>\n",
              "      <th>GC</th>\n",
              "    </tr>\n",
              "  </thead>\n",
              "  <tbody>\n",
              "    <tr>\n",
              "      <th>0</th>\n",
              "      <td>chr1</td>\n",
              "      <td>0</td>\n",
              "      <td>20000</td>\n",
              "      <td>NaN</td>\n",
              "      <td>0.59230</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>1</th>\n",
              "      <td>chr1</td>\n",
              "      <td>20000</td>\n",
              "      <td>40000</td>\n",
              "      <td>NaN</td>\n",
              "      <td>0.49430</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>2</th>\n",
              "      <td>chr1</td>\n",
              "      <td>40000</td>\n",
              "      <td>60000</td>\n",
              "      <td>NaN</td>\n",
              "      <td>0.36075</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>3</th>\n",
              "      <td>chr1</td>\n",
              "      <td>60000</td>\n",
              "      <td>80000</td>\n",
              "      <td>NaN</td>\n",
              "      <td>0.35660</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>4</th>\n",
              "      <td>chr1</td>\n",
              "      <td>80000</td>\n",
              "      <td>100000</td>\n",
              "      <td>NaN</td>\n",
              "      <td>0.40255</td>\n",
              "    </tr>\n",
              "  </tbody>\n",
              "</table>\n",
              "</div>"
            ],
            "text/plain": [
              "  chrom  start     end  weight       GC\n",
              "0  chr1      0   20000     NaN  0.59230\n",
              "1  chr1  20000   40000     NaN  0.49430\n",
              "2  chr1  40000   60000     NaN  0.36075\n",
              "3  chr1  60000   80000     NaN  0.35660\n",
              "4  chr1  80000  100000     NaN  0.40255"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 51
        }
      ]
    },
    {
      "metadata": {
        "id": "-3dbRJVHqKXY",
        "colab_type": "code",
        "colab": {}
      },
      "cell_type": "code",
      "source": [
        "lam, vec = cooltools.eigdecomp.cooler_cis_eig(\n",
        "    clr, \n",
        "    bins, \n",
        "    regions=[x[:3] for x in chromarms.itertuples(index=False)], \n",
        "    n_eigs=3, \n",
        "    phasing_track_col='GC', \n",
        "    balance='weight', \n",
        "    ignore_diags=None, \n",
        "    clip_percentile=99.9, \n",
        "    sort_metric=None\n",
        ")"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "metadata": {
        "id": "t55HR7RbsjED",
        "colab_type": "code",
        "outputId": "cc44786e-7748-4dfb-f753-3943929f35fd",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 204
        }
      },
      "cell_type": "code",
      "source": [
        "vec.head()"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/html": [
              "<div>\n",
              "<style scoped>\n",
              "    .dataframe tbody tr th:only-of-type {\n",
              "        vertical-align: middle;\n",
              "    }\n",
              "\n",
              "    .dataframe tbody tr th {\n",
              "        vertical-align: top;\n",
              "    }\n",
              "\n",
              "    .dataframe thead th {\n",
              "        text-align: right;\n",
              "    }\n",
              "</style>\n",
              "<table border=\"1\" class=\"dataframe\">\n",
              "  <thead>\n",
              "    <tr style=\"text-align: right;\">\n",
              "      <th></th>\n",
              "      <th>chrom</th>\n",
              "      <th>start</th>\n",
              "      <th>end</th>\n",
              "      <th>weight</th>\n",
              "      <th>GC</th>\n",
              "      <th>E1</th>\n",
              "      <th>E2</th>\n",
              "      <th>E3</th>\n",
              "    </tr>\n",
              "  </thead>\n",
              "  <tbody>\n",
              "    <tr>\n",
              "      <th>0</th>\n",
              "      <td>chr1</td>\n",
              "      <td>0</td>\n",
              "      <td>20000</td>\n",
              "      <td>NaN</td>\n",
              "      <td>0.59230</td>\n",
              "      <td>NaN</td>\n",
              "      <td>NaN</td>\n",
              "      <td>NaN</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>1</th>\n",
              "      <td>chr1</td>\n",
              "      <td>20000</td>\n",
              "      <td>40000</td>\n",
              "      <td>NaN</td>\n",
              "      <td>0.49430</td>\n",
              "      <td>NaN</td>\n",
              "      <td>NaN</td>\n",
              "      <td>NaN</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>2</th>\n",
              "      <td>chr1</td>\n",
              "      <td>40000</td>\n",
              "      <td>60000</td>\n",
              "      <td>NaN</td>\n",
              "      <td>0.36075</td>\n",
              "      <td>NaN</td>\n",
              "      <td>NaN</td>\n",
              "      <td>NaN</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>3</th>\n",
              "      <td>chr1</td>\n",
              "      <td>60000</td>\n",
              "      <td>80000</td>\n",
              "      <td>NaN</td>\n",
              "      <td>0.35660</td>\n",
              "      <td>NaN</td>\n",
              "      <td>NaN</td>\n",
              "      <td>NaN</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>4</th>\n",
              "      <td>chr1</td>\n",
              "      <td>80000</td>\n",
              "      <td>100000</td>\n",
              "      <td>NaN</td>\n",
              "      <td>0.40255</td>\n",
              "      <td>NaN</td>\n",
              "      <td>NaN</td>\n",
              "      <td>NaN</td>\n",
              "    </tr>\n",
              "  </tbody>\n",
              "</table>\n",
              "</div>"
            ],
            "text/plain": [
              "  chrom  start     end  weight       GC  E1  E2  E3\n",
              "0  chr1      0   20000     NaN  0.59230 NaN NaN NaN\n",
              "1  chr1  20000   40000     NaN  0.49430 NaN NaN NaN\n",
              "2  chr1  40000   60000     NaN  0.36075 NaN NaN NaN\n",
              "3  chr1  60000   80000     NaN  0.35660 NaN NaN NaN\n",
              "4  chr1  80000  100000     NaN  0.40255 NaN NaN NaN"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 55
        }
      ]
    },
    {
      "metadata": {
        "id": "-RoHvM5_2P4t",
        "colab_type": "code",
        "outputId": "2eefb16d-d073-4dd2-b258-91c25f7c1801",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 235
        }
      },
      "cell_type": "code",
      "source": [
        "lam.head()"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/html": [
              "<div>\n",
              "<style scoped>\n",
              "    .dataframe tbody tr th:only-of-type {\n",
              "        vertical-align: middle;\n",
              "    }\n",
              "\n",
              "    .dataframe tbody tr th {\n",
              "        vertical-align: top;\n",
              "    }\n",
              "\n",
              "    .dataframe thead th {\n",
              "        text-align: right;\n",
              "    }\n",
              "</style>\n",
              "<table border=\"1\" class=\"dataframe\">\n",
              "  <thead>\n",
              "    <tr style=\"text-align: right;\">\n",
              "      <th></th>\n",
              "      <th>eigval1</th>\n",
              "      <th>eigval2</th>\n",
              "      <th>eigval3</th>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>region</th>\n",
              "      <th></th>\n",
              "      <th></th>\n",
              "      <th></th>\n",
              "    </tr>\n",
              "  </thead>\n",
              "  <tbody>\n",
              "    <tr>\n",
              "      <th>chr1:0-123035434</th>\n",
              "      <td>1715.509511</td>\n",
              "      <td>1119.805776</td>\n",
              "      <td>-890.544007</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>chr1:123035434-249250621</th>\n",
              "      <td>1142.952018</td>\n",
              "      <td>1026.974430</td>\n",
              "      <td>-868.326597</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>chr2:0-93826171</th>\n",
              "      <td>1382.297385</td>\n",
              "      <td>770.204333</td>\n",
              "      <td>-690.263584</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>chr2:93826171-243199373</th>\n",
              "      <td>1845.107407</td>\n",
              "      <td>1106.526004</td>\n",
              "      <td>-1104.508453</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>chr3:0-92004854</th>\n",
              "      <td>1392.531753</td>\n",
              "      <td>807.529687</td>\n",
              "      <td>-800.195451</td>\n",
              "    </tr>\n",
              "  </tbody>\n",
              "</table>\n",
              "</div>"
            ],
            "text/plain": [
              "                              eigval1      eigval2      eigval3\n",
              "region                                                         \n",
              "chr1:0-123035434          1715.509511  1119.805776  -890.544007\n",
              "chr1:123035434-249250621  1142.952018  1026.974430  -868.326597\n",
              "chr2:0-93826171           1382.297385   770.204333  -690.263584\n",
              "chr2:93826171-243199373   1845.107407  1106.526004 -1104.508453\n",
              "chr3:0-92004854           1392.531753   807.529687  -800.195451"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 92
        }
      ]
    },
    {
      "metadata": {
        "id": "91NkmcQGB9wp",
        "colab_type": "text"
      },
      "cell_type": "markdown",
      "source": [
        "### Write and download file to local drive"
      ]
    },
    {
      "metadata": {
        "id": "SyJACHsq3Ily",
        "colab_type": "code",
        "outputId": "af50e27d-8e2e-49de-d187-6848b0c0d15a",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 34
        }
      },
      "cell_type": "code",
      "source": [
        "#write dataframe to file\n",
        "bioframe.to_bigwig(vec, hg19, 'eigen_byarm.bw', 'E1')"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "bedGraphToBigWig /tmp/tmpsprvbqnp.bg /tmp/tmpo0112kb7.chrom.sizes eigen_byarm.bw\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "metadata": {
        "id": "nswV0GP78BHB",
        "colab_type": "code",
        "colab": {}
      },
      "cell_type": "code",
      "source": [
        "from google.colab import files\n",
        "\n",
        "files.download('eigen_byarm.bw')"
      ],
      "execution_count": 0,
      "outputs": []
    }
  ]
}