Datasets
Models
Downloads: 1
1022 lines (1021 with data), 63.1 kB
{
"cells": [
{
"cell_type": "code",
"execution_count": 63,
"metadata": {},
"outputs": [],
"source": [
"%matplotlib inline\n",
"import matplotlib as mpl\n",
"mpl.rcParams['figure.dpi'] = 250\n",
"\n",
"import matplotlib.pyplot as plt\n",
"\n",
"import pysam\n",
"import singlecellmultiomics.molecule\n",
"import singlecellmultiomics.fragment\n",
"import pysamiterators\n",
"import pandas as pd\n",
"\n",
"nla_test_bam_path = '../data/mini_nla_test.bam'"
]
},
{
"cell_type": "code",
"execution_count": 64,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"ST-E00285:221:HKM7VCCXY:6:2211:11485:25534\t1171\t0\t164834715\t60\t151M\t0\t164834656\t151\tGAGTTATGAAATCCCTTGCTAACTTTCTCTTCTCTGGGGAAAAGAGTCTGAGTTCCCTCAGCCTTTCTTCCTAAGACCTGTGGTTGGCTATAAATTGCATTGGTTGCTCTTCTCTAATGACCTCTGGAGATTGGCAACACCAAACAAAGAC\tarray('B', [12, 41, 41, 32, 41, 37, 37, 37, 41, 41, 41, 37, 37, 41, 37, 32, 27, 37, 41, 41, 41, 32, 37, 32, 41, 41, 41, 41, 41, 37, 37, 37, 41, 41, 41, 37, 32, 37, 37, 41, 37, 41, 41, 41, 41, 41, 41, 41, 41, 41, 37, 37, 27, 41, 41, 27, 41, 32, 41, 41, 41, 37, 41, 41, 41, 41, 37, 37, 37, 41, 41, 41, 41, 37, 27, 41, 41, 41, 41, 41, 41, 32, 37, 22, 32, 41, 32, 27, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 37, 41, 41, 41, 37, 37, 41, 41, 41, 41, 37, 41, 41, 41, 41, 37, 41, 41, 41, 41, 27, 37, 41, 41, 41, 41, 37, 41, 41, 41, 37, 41, 37, 37, 41, 37, 37, 41, 41, 37, 41, 41, 41, 37, 41, 37, 37, 37, 32, 32])\t[('NM', 0), ('MD', '151'), ('MC', '140M'), ('AS', 151), ('XS', 20), ('Is', 'ST-E00285'), ('RN', '221'), ('Fc', 'HKM7VCCXY'), ('La', '6'), ('Ti', '2211'), ('CX', '11485'), ('CY', '25534'), ('Fi', 'N'), ('CN', '0'), ('aA', 'TATAATAT'), ('LY', 'APKS1P25-NLAP1L4'), ('RX', 'TTA'), ('RQ', 'AAF'), ('BI', '90'), ('bc', 'AGCAAGCA'), ('BC', 'AGCAAGCA'), ('QT', 'FFFJJJJJ'), ('MX', 'NLAIII96C8U3'), ('MI', 'AGCAAGCATTATATAATAT'), ('QM', 'FFFJJJJJAAF////////'), ('SM', 'APKS1P25-NLAP1L4_90'), ('rS', 'GTCTTT'), ('rP', 164834866), ('RC', 2), ('DS', 164834656), ('DT', 'NLA'), ('RS', '-'), ('RZ', 'CATG'), ('RG', 'HKM7VCCXY.6.APKS1P25-NLAP1L4_90')]\n",
"ST-E00285:221:HKM7VCCXY:5:1211:30523:24796\t147\t0\t164834716\t60\t151M\t0\t164834656\t151\tAGTTATGAAATCCCTTGCTAACTTTCTCTTCTCTGGGGAAAAGAGTCTGAGTTCCCTCAGCCTTTCTTCCTAAGACCTGTGGTTGGCTATAAATTGCATTGGTTGCTCTTCTCTAATGACCTCTGGAGATTGGCAACACCAAACAAAAACA\tarray('B', [32, 41, 32, 27, 41, 41, 37, 37, 22, 37, 27, 32, 37, 37, 32, 41, 41, 41, 41, 41, 41, 41, 41, 27, 37, 41, 41, 37, 32, 32, 37, 27, 27, 41, 41, 41, 41, 41, 37, 37, 27, 37, 32, 37, 37, 37, 32, 37, 37, 27, 32, 37, 12, 32, 37, 22, 41, 41, 41, 32, 22, 22, 37, 22, 32, 37, 41, 32, 27, 32, 22, 32, 12, 41, 41, 37, 37, 41, 41, 41, 41, 41, 41, 37, 22, 12, 22, 12, 32, 41, 41, 32, 37, 32, 27, 12, 12, 37, 37, 22, 32, 12, 37, 37, 12, 32, 27, 12, 22, 32, 27, 41, 37, 37, 37, 32, 41, 41, 41, 41, 27, 37, 41, 32, 37, 27, 12, 41, 12, 37, 22, 37, 32, 22, 32, 12, 37, 32, 12, 41, 41, 32, 12, 32, 27, 37, 12, 32, 32, 12, 27])\t[('NM', 1), ('MD', '147G3'), ('MC', '140M'), ('AS', 147), ('XS', 20), ('Is', 'ST-E00285'), ('RN', '221'), ('Fc', 'HKM7VCCXY'), ('La', '5'), ('Ti', '1211'), ('CX', '30523'), ('CY', '24796'), ('Fi', 'N'), ('CN', '0'), ('aA', 'TTAGGCAT'), ('LY', 'APKS1-P10-1-3'), ('RX', 'ACG'), ('RQ', 'AAF'), ('BI', '29'), ('bc', 'CGTTGATC'), ('BC', 'CGTTGATC'), ('QT', 'FFJJJJFA'), ('MX', 'NLAIII96C8U3'), ('MI', 'CGTTGATCACGTTAGGCAT'), ('QM', 'FFJJJJFAAAF////////'), ('SM', 'APKS1-P10-1-3_29'), ('rS', 'TGTTTT'), ('rP', 164834867), ('RC', 1), ('DS', 164834656), ('DT', 'NLA'), ('RS', '-'), ('RZ', 'CATG'), ('RG', 'HKM7VCCXY.5.APKS1-P10-1-3_29')]\n"
]
}
],
"source": [
"# iterate over the reads in nla_test_bam_path using pysam ...\n",
"with pysam.AlignmentFile(nla_test_bam_path) as alignments:\n",
" for i,read in enumerate( alignments ):\n",
" print(str(read))\n",
" if i>=1:\n",
" break"
]
},
{
"cell_type": "code",
"execution_count": 65,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"ST-E00285:221:HKM7VCCXY:5:2221:19918:29490\t83\t0\t164834729\t60\t140M\t0\t164834726\t140\tCTTGCTAACTTTCTCTTCTCTGGGGAAAAGAGTCTGAGTTCCCTCAGCCTTTCTTCCTAAGACCTGTGGTTGGCTATAAATTGCATTGGTTGCTCTTCTCTAATGACCTCTGGAGATTGGCAACACCAAACAAAGACATG\tarray('B', [37, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 37, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 37, 37, 32, 32, 37, 37, 37, 41, 41, 37, 41, 41, 41, 41, 37, 41, 32, 32, 32, 22, 37, 12, 41, 37, 41, 41, 27, 41, 37, 37, 27, 37, 37, 37, 27, 37, 41, 37, 37, 37, 41, 41, 37, 37, 32, 41, 32, 32, 37, 32])\t[('NM', 0), ('MD', '140'), ('MC', '143M8S'), ('AS', 140), ('XS', 20), ('Is', 'ST-E00285'), ('RN', '221'), ('Fc', 'HKM7VCCXY'), ('La', '5'), ('Ti', '2221'), ('CX', '19918'), ('CY', '29490'), ('Fi', 'N'), ('CN', '0'), ('aA', 'TTAGGCAT'), ('LY', 'APKS1-P10-1-3'), ('RX', 'GAC'), ('RQ', 'AA<'), ('BI', '75'), ('bc', 'CAAGTTCC'), ('BC', 'CAAGTTCC'), ('QT', 'AAFAAFFJ'), ('MX', 'NLAIII96C8U3'), ('MI', 'CAAGTTCCGACTTAGGCAT'), ('QM', 'AAFAAFFJAA<////////'), ('SM', 'APKS1-P10-1-3_75'), ('rS', 'TCCCTT'), ('rP', 164834726), ('RC', 1), ('DS', 164834865), ('DT', 'NLA'), ('RS', '+'), ('RZ', 'CATG'), ('RG', 'HKM7VCCXY.5.APKS1-P10-1-3_75')] ST-E00285:221:HKM7VCCXY:5:2221:19918:29490\t163\t0\t164834726\t60\t143M8S\t0\t164834729\t143\tTCCCTTGCTAACTTTCTCTTCTCTGGGGAAAAGAGTCTGACTTCCCTCAGCCTTTCTTCCTAAGACCTGTGGTTGGCTATAAATTGCATTGGTTGCTCTTCTCTAATGACCTCTGGAGATTGGCAACACCAAACAAAGACATGGGAACTTG\tarray('B', [32, 32, 27, 27, 12, 27, 32, 37, 41, 27, 12, 22, 12, 22, 27, 37, 37, 12, 22, 37, 12, 12, 37, 32, 37, 12, 22, 12, 12, 32, 41, 12, 37, 12, 37, 41, 27, 22, 37, 27, 12, 27, 12, 27, 27, 32, 22, 32, 37, 22, 32, 37, 41, 12, 37, 37, 41, 37, 32, 27, 37, 41, 41, 37, 41, 41, 41, 41, 37, 37, 41, 41, 37, 37, 41, 37, 41, 27, 37, 22, 41, 37, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 32, 41, 41, 41, 32, 41, 41, 41, 41, 41, 41, 41, 41, 32, 32, 41, 41, 37, 41, 37, 37, 37, 22, 27, 37, 41, 41, 37, 37, 37, 41, 41, 27, 37, 37, 41, 41, 41, 41, 27, 37, 37, 37, 37])\t[('NM', 1), ('MD', '40G102'), ('MC', '140M'), ('AS', 138), ('XS', 20), ('Is', 'ST-E00285'), ('RN', '221'), ('Fc', 'HKM7VCCXY'), ('La', '5'), ('Ti', '2221'), ('CX', '19918'), ('CY', '29490'), ('Fi', 'N'), ('CN', '0'), ('aA', 'TTAGGCAT'), ('LY', 'APKS1-P10-1-3'), ('RX', 'GAC'), ('RQ', 'AA<'), ('BI', '75'), ('bc', 'CAAGTTCC'), ('BC', 'CAAGTTCC'), ('QT', 'AAFAAFFJ'), ('MX', 'NLAIII96C8U3'), ('MI', 'CAAGTTCCGACTTAGGCAT'), ('QM', 'AAFAAFFJAA<////////'), ('SM', 'APKS1-P10-1-3_75'), ('rS', 'TCCCTT'), ('rP', 164834726), ('RC', 1), ('DS', 164834865), ('DT', 'NLA'), ('RS', '+'), ('RZ', 'CATG'), ('RG', 'HKM7VCCXY.5.APKS1-P10-1-3_75')]\n",
"ST-E00285:221:HKM7VCCXY:6:1214:26920:26607\t83\t0\t164834729\t60\t140M\t0\t164834729\t140\tCTTGCTAACTTTCTCTTCTCTGGGGAAAAGAGTCTGAGTTCCCTCAGCCTTTCTTCCTAAGACCTGTGGTTGGCTATAAATTGCATTGGTTGCTCTTCTCTAATGACCTCTGGAGATTGGCAACACCAAACAAAGACATG\tarray('B', [37, 37, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 27, 41, 41, 27, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 37, 41, 41, 41, 41, 41, 41, 41, 41, 32, 41, 41, 41, 37, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 37, 41, 41, 41, 41, 41, 32, 37, 37, 41, 41, 41, 41, 41, 41, 41, 41, 37, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 37, 41, 41, 37, 41, 41, 41, 37, 37, 41, 37, 22, 32, 41, 41, 41, 41, 37, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41])\t[('NM', 0), ('MD', '140'), ('MC', '140M11S'), ('AS', 140), ('XS', 20), ('Is', 'ST-E00285'), ('RN', '221'), ('Fc', 'HKM7VCCXY'), ('La', '6'), ('Ti', '1214'), ('CX', '26920'), ('CY', '26607'), ('Fi', 'N'), ('CN', '0'), ('aA', 'CTCAGAAT'), ('LY', 'APKS1P25-NLAP1L1'), ('RX', 'ACC'), ('RQ', 'AAA'), ('BI', '61'), ('bc', 'CCTTCACA'), ('BC', 'CCTTCACA'), ('QT', 'FFJFJJJJ'), ('MX', 'NLAIII96C8U3'), ('MI', 'CCTTCACAACCCTCAGAAT'), ('QM', 'FFJFJJJJAAA////////'), ('SM', 'APKS1P25-NLAP1L1_61'), ('rS', 'CTTGCT'), ('rP', 164834729), ('RC', 1), ('DS', 164834865), ('DT', 'NLA'), ('RS', '+'), ('RZ', 'CATG'), ('RG', 'HKM7VCCXY.6.APKS1P25-NLAP1L1_61')] ST-E00285:221:HKM7VCCXY:6:1214:26920:26607\t163\t0\t164834729\t60\t140M11S\t0\t164834729\t140\tCTTGCTAACTTTCTCTTCTCTGGGGAAAAGAGTCTGAGTTCCCTCAGCCTTTCTTCCTAAGACCTGTGGTTGGCTATAAATTGCATTGGTTGCTCTTCTCTAATGACCTCTGGAGATTGGCAACACCAAACAAAGACATGTGTGAAGGGGT\tarray('B', [32, 32, 37, 37, 37, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 37, 37, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 37, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 37, 37, 41, 41, 41, 41, 37, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 27, 37, 41, 41, 37, 41, 41, 41, 41, 41, 41, 22, 37, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 37, 41, 37, 41, 37, 41, 41, 37, 37, 41, 41, 41, 12])\t[('NM', 0), ('MD', '140'), ('MC', '140M'), ('AS', 140), ('XS', 20), ('Is', 'ST-E00285'), ('RN', '221'), ('Fc', 'HKM7VCCXY'), ('La', '6'), ('Ti', '1214'), ('CX', '26920'), ('CY', '26607'), ('Fi', 'N'), ('CN', '0'), ('aA', 'CTCAGAAT'), ('LY', 'APKS1P25-NLAP1L1'), ('RX', 'ACC'), ('RQ', 'AAA'), ('BI', '61'), ('bc', 'CCTTCACA'), ('BC', 'CCTTCACA'), ('QT', 'FFJFJJJJ'), ('MX', 'NLAIII96C8U3'), ('MI', 'CCTTCACAACCCTCAGAAT'), ('QM', 'FFJFJJJJAAA////////'), ('SM', 'APKS1P25-NLAP1L1_61'), ('rS', 'CTTGCT'), ('rP', 164834729), ('RC', 1), ('DS', 164834865), ('DT', 'NLA'), ('RS', '+'), ('RZ', 'CATG'), ('RG', 'HKM7VCCXY.6.APKS1P25-NLAP1L1_61')]\n"
]
}
],
"source": [
"# Iterate over the mate pairs, keeping R1 and R2 together\n",
"with pysam.AlignmentFile(nla_test_bam_path) as alignments:\n",
" for i,(R1,R2) in enumerate( pysamiterators.MatePairIterator(alignments) ):\n",
" print(R1, R2)\n",
" if i>=1:\n",
" break"
]
},
{
"cell_type": "code",
"execution_count": 66,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Molecule\n",
" with 4 assinged fragments\n",
" Allele :No allele assigned\n",
" Fragment:\n",
" sample:APKS1P25-NLAP2L1_30\n",
" umi:CCG\n",
" span:chr1 164834731-164834869\n",
" strand:-\n",
" has R1: no\n",
" has R2: yes\n",
" \n",
"\t Fragment:\n",
" sample:APKS1P25-NLAP2L1_30\n",
" umi:CCG\n",
" span:chr1 164834731-164834869\n",
" strand:-\n",
" has R1: no\n",
" has R2: yes\n",
" \n",
"\t Fragment:\n",
" sample:APKS1P25-NLAP2L1_30\n",
" umi:CCG\n",
" span:chr1 164834731-164834869\n",
" strand:-\n",
" has R1: no\n",
" has R2: yes\n",
" \n",
"\t Fragment:\n",
" sample:APKS1P25-NLAP2L1_30\n",
" umi:CCG\n",
" span:chr1 164834741-164834869\n",
" strand:-\n",
" has R1: no\n",
" has R2: yes\n",
" \n"
]
}
],
"source": [
"# Iterate over molecules\n",
"with pysam.AlignmentFile(nla_test_bam_path) as alignments:\n",
" for i,molecule in enumerate( singlecellmultiomics.molecule.MoleculeIterator(alignments) ):\n",
" print(molecule)\n",
" break"
]
},
{
"cell_type": "code",
"execution_count": 67,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"Fragment:\n",
" sample:APKS1P25-NLAP2L1_30\n",
" umi:CCG\n",
" span:chr1 164834731-164834869\n",
" strand:-\n",
" has R1: no\n",
" has R2: yes\n",
" "
]
},
"execution_count": 67,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Obtain the first fragment assigned to the molecule\n",
"molecule[0]"
]
},
{
"cell_type": "code",
"execution_count": 68,
"metadata": {},
"outputs": [],
"source": [
"# obtain read 1 from the first fragment\n",
"molecule[0][0] "
]
},
{
"cell_type": "code",
"execution_count": 69,
"metadata": {},
"outputs": [],
"source": [
"molecule[0].get_R1()"
]
},
{
"cell_type": "code",
"execution_count": 70,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"'-'"
]
},
"execution_count": 70,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Obtain the strand as string (-: reverse, +:forward)\n",
"molecule.get_strand_repr()"
]
},
{
"cell_type": "code",
"execution_count": 71,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"4"
]
},
"execution_count": 71,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Obtain the amount of fragments associated to the molecule\n",
"len(molecule)"
]
},
{
"cell_type": "code",
"execution_count": 72,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"'APKS1P25-NLAP2L1_30'"
]
},
"execution_count": 72,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Obtain the sample name\n",
"molecule.sample"
]
},
{
"cell_type": "code",
"execution_count": 73,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"NlaIIIMolecule\n",
" with 4 assinged fragments\n",
" Allele :No allele assigned\n",
" Fragment:\n",
" sample:APKS3-P19-1-1_115\n",
" umi:GAC\n",
" span:chr1 164834869-164835153\n",
" strand:-\n",
" has R1: yes\n",
" has R2: yes\n",
" DS:164834865\n",
" \tRS:1\n",
" \tRZ:CATG\n",
" \tRestriction site:('chr1', 164834865)\n",
"\t Fragment:\n",
" sample:APKS3-P19-1-1_115\n",
" umi:GCC\n",
" span:chr1 164834869-164835270\n",
" strand:-\n",
" has R1: yes\n",
" has R2: yes\n",
" DS:164834865\n",
" \tRS:1\n",
" \tRZ:CATG\n",
" \tRestriction site:('chr1', 164834865)\n",
"\t Fragment:\n",
" sample:APKS3-P19-1-1_115\n",
" umi:GAC\n",
" span:chr1 164834869-164835270\n",
" strand:-\n",
" has R1: yes\n",
" has R2: yes\n",
" DS:164834865\n",
" \tRS:1\n",
" \tRZ:CATG\n",
" \tRestriction site:('chr1', 164834865)\n",
"\t Fragment:\n",
" sample:APKS3-P19-1-1_115\n",
" umi:GAC\n",
" span:chr1 164834869-164835270\n",
" strand:-\n",
" has R1: yes\n",
" has R2: yes\n",
" DS:164834865\n",
" \tRS:1\n",
" \tRZ:CATG\n",
" \tRestriction site:('chr1', 164834865)"
]
},
"execution_count": 73,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Find a molecule with 4 fragments\n",
"molecules_seen = [] # store all molecules, used in next examples\n",
"with pysam.AlignmentFile(nla_test_bam_path) as alignments:\n",
" for i,molecule in enumerate(\n",
" singlecellmultiomics.molecule.MoleculeIterator(alignments,\n",
" fragment_class_args={\n",
" 'umi_hamming_distance':1\n",
" },\n",
" moleculeClass=singlecellmultiomics.molecule.NlaIIIMolecule,\n",
" fragmentClass=singlecellmultiomics.fragment.NLAIIIFragment\n",
"\n",
" )):\n",
" molecules_seen.append(molecule)\n",
" if len(molecule)==4 and i>0:\n",
" break\n",
"molecule"
]
},
{
"cell_type": "code",
"execution_count": 74,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"Counter({'GAC': 3, 'GCC': 1})"
]
},
"execution_count": 74,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Obtain associated unique molecular identifiers\n",
"molecule.umi_counter"
]
},
{
"cell_type": "code",
"execution_count": 75,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"ST-E00285:221:HKM7VCCXY:6:1224:9120:30439\n",
"NS500413:404:HJ3KHBGX5:2:22312:23759:6537\n",
"NS500413:404:HJ3KHBGX5:4:13507:22698:18344\n",
"NS500413:404:HJ3KHBGX5:3:23502:18585:5263\n"
]
}
],
"source": [
"# Iterate over all fragments in the molecule, obtain their R1 and print the read name:\n",
"for fragment in molecule:\n",
" print(fragment.get_R1().query_name )"
]
},
{
"cell_type": "code",
"execution_count": 76,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"164835159 TGCAGT [Fragment:\n",
" sample:APKS3-P19-1-1_115\n",
" umi:GAC\n",
" span:chr1 164834869-164835153\n",
" strand:-\n",
" has R1: yes\n",
" has R2: yes\n",
" DS:164834865\n",
"\tRS:1\n",
"\tRZ:CATG\n",
"\tRestriction site:('chr1', 164834865)]\n",
"164835276 CAGTGT [Fragment:\n",
" sample:APKS3-P19-1-1_115\n",
" umi:GCC\n",
" span:chr1 164834869-164835270\n",
" strand:-\n",
" has R1: yes\n",
" has R2: yes\n",
" DS:164834865\n",
"\tRS:1\n",
"\tRZ:CATG\n",
"\tRestriction site:('chr1', 164834865)]\n",
"164835276 CCGTGT [Fragment:\n",
" sample:APKS3-P19-1-1_115\n",
" umi:GAC\n",
" span:chr1 164834869-164835270\n",
" strand:-\n",
" has R1: yes\n",
" has R2: yes\n",
" DS:164834865\n",
"\tRS:1\n",
"\tRZ:CATG\n",
"\tRestriction site:('chr1', 164834865), Fragment:\n",
" sample:APKS3-P19-1-1_115\n",
" umi:GAC\n",
" span:chr1 164834869-164835270\n",
" strand:-\n",
" has R1: yes\n",
" has R2: yes\n",
" DS:164834865\n",
"\tRS:1\n",
"\tRZ:CATG\n",
"\tRestriction site:('chr1', 164834865)]\n"
]
}
],
"source": [
"# Iterate over all reverse transcription reactions:\n",
"for (reverse_primer_start, reverse_primer_sequence), associated_fragments in molecule.get_rt_reactions().items():\n",
" print(reverse_primer_start,reverse_primer_sequence,associated_fragments)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Equivalence testing"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Comparing fragments"
]
},
{
"cell_type": "code",
"execution_count": 77,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"True"
]
},
"execution_count": 77,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# compare two fragments: (check if they should belong to the same molecule)\n",
"fragment_A = molecule[0]\n",
"fragment_B = molecule[1]\n",
"fragment_A == fragment_B"
]
},
{
"cell_type": "code",
"execution_count": 78,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"Fragment:\n",
" sample:APKS1P25-NLAP2L2_57\n",
" umi:CCG\n",
" span:chr1 164834728-164834868\n",
" strand:+\n",
" has R1: yes\n",
" has R2: no\n",
" DS:164834865\n",
"\tRS:0\n",
"\tRZ:CAT\n",
"\tRestriction site:('chr1', 164834865)"
]
},
"execution_count": 78,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Obtain a fragment not belonging to the molecule \n",
"fragment_C = molecules_seen[0][0]\n",
"fragment_C"
]
},
{
"cell_type": "code",
"execution_count": 79,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"False"
]
},
"execution_count": 79,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"fragment_C == fragment_A"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Comparing fragment to molecule"
]
},
{
"cell_type": "code",
"execution_count": 80,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"True"
]
},
"execution_count": 80,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Fragment A belongs to molecule, this comparison results in True\n",
"fragment_A == molecule"
]
},
{
"cell_type": "code",
"execution_count": 81,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"False"
]
},
"execution_count": 81,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Fragment C does not belong to molecule\n",
"fragment_C == molecule "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Consensus sequence"
]
},
{
"cell_type": "code",
"execution_count": 82,
"metadata": {},
"outputs": [
{
"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 tr th {\n",
" text-align: left;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr>\n",
" <th></th>\n",
" <th colspan=\"21\" halign=\"left\">chr1</th>\n",
" </tr>\n",
" <tr>\n",
" <th></th>\n",
" <th>164834869</th>\n",
" <th>164834870</th>\n",
" <th>164834871</th>\n",
" <th>164834872</th>\n",
" <th>164834873</th>\n",
" <th>164834874</th>\n",
" <th>164834875</th>\n",
" <th>164834876</th>\n",
" <th>164834877</th>\n",
" <th>164834878</th>\n",
" <th>...</th>\n",
" <th>164835262</th>\n",
" <th>164835263</th>\n",
" <th>164835264</th>\n",
" <th>164835265</th>\n",
" <th>164835266</th>\n",
" <th>164835267</th>\n",
" <th>164835268</th>\n",
" <th>164835269</th>\n",
" <th>164835270</th>\n",
" <th>164835238</th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <td>A</td>\n",
" <td>4.0</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>4.0</td>\n",
" <td>NaN</td>\n",
" <td>4.0</td>\n",
" <td>NaN</td>\n",
" <td>4.0</td>\n",
" <td>NaN</td>\n",
" <td>...</td>\n",
" <td>NaN</td>\n",
" <td>3.0</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>3.0</td>\n",
" <td>3.0</td>\n",
" <td>NaN</td>\n",
" <td>3.0</td>\n",
" <td>3.0</td>\n",
" <td>NaN</td>\n",
" </tr>\n",
" <tr>\n",
" <td>G</td>\n",
" <td>NaN</td>\n",
" <td>4.0</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>4.0</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>...</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>1.0</td>\n",
" </tr>\n",
" <tr>\n",
" <td>T</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>4.0</td>\n",
" <td>4.0</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>4.0</td>\n",
" <td>NaN</td>\n",
" <td>4.0</td>\n",
" <td>...</td>\n",
" <td>3.0</td>\n",
" <td>NaN</td>\n",
" <td>3.0</td>\n",
" <td>3.0</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>3.0</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" </tr>\n",
" <tr>\n",
" <td>C</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>...</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"<p>4 rows × 353 columns</p>\n",
"</div>"
],
"text/plain": [
" chr1 \\\n",
" 164834869 164834870 164834871 164834872 164834873 164834874 164834875 \n",
"A 4.0 NaN NaN NaN 4.0 NaN 4.0 \n",
"G NaN 4.0 NaN NaN NaN 4.0 NaN \n",
"T NaN NaN 4.0 4.0 NaN NaN NaN \n",
"C NaN NaN NaN NaN NaN NaN NaN \n",
"\n",
" ... \\\n",
" 164834876 164834877 164834878 ... 164835262 164835263 164835264 164835265 \n",
"A NaN 4.0 NaN ... NaN 3.0 NaN NaN \n",
"G NaN NaN NaN ... NaN NaN NaN NaN \n",
"T 4.0 NaN 4.0 ... 3.0 NaN 3.0 3.0 \n",
"C NaN NaN NaN ... NaN NaN NaN NaN \n",
"\n",
" \n",
" 164835266 164835267 164835268 164835269 164835270 164835238 \n",
"A 3.0 3.0 NaN 3.0 3.0 NaN \n",
"G NaN NaN NaN NaN NaN 1.0 \n",
"T NaN NaN 3.0 NaN NaN NaN \n",
"C NaN NaN NaN NaN NaN NaN \n",
"\n",
"[4 rows x 353 columns]"
]
},
"execution_count": 82,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Obtain the aligned base frequencies in a pandas dataframe\n",
"pd.DataFrame( molecule.get_base_observation_dict() )"
]
},
{
"cell_type": "code",
"execution_count": 83,
"metadata": {},
"outputs": [
{
"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 tr th {\n",
" text-align: left;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr>\n",
" <th></th>\n",
" <th colspan=\"21\" halign=\"left\">chr1</th>\n",
" </tr>\n",
" <tr>\n",
" <th></th>\n",
" <th>164834869</th>\n",
" <th>164834870</th>\n",
" <th>164834871</th>\n",
" <th>164834872</th>\n",
" <th>164834873</th>\n",
" <th>164834874</th>\n",
" <th>164834875</th>\n",
" <th>164834876</th>\n",
" <th>164834877</th>\n",
" <th>164834878</th>\n",
" <th>...</th>\n",
" <th>164835261</th>\n",
" <th>164835262</th>\n",
" <th>164835263</th>\n",
" <th>164835264</th>\n",
" <th>164835265</th>\n",
" <th>164835266</th>\n",
" <th>164835267</th>\n",
" <th>164835268</th>\n",
" <th>164835269</th>\n",
" <th>164835270</th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <td>base</td>\n",
" <td>A</td>\n",
" <td>G</td>\n",
" <td>T</td>\n",
" <td>T</td>\n",
" <td>A</td>\n",
" <td>G</td>\n",
" <td>A</td>\n",
" <td>T</td>\n",
" <td>A</td>\n",
" <td>T</td>\n",
" <td>...</td>\n",
" <td>T</td>\n",
" <td>T</td>\n",
" <td>A</td>\n",
" <td>T</td>\n",
" <td>T</td>\n",
" <td>A</td>\n",
" <td>A</td>\n",
" <td>T</td>\n",
" <td>A</td>\n",
" <td>A</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"<p>1 rows × 353 columns</p>\n",
"</div>"
],
"text/plain": [
" chr1 \\\n",
" 164834869 164834870 164834871 164834872 164834873 164834874 164834875 \n",
"base A G T T A G A \n",
"\n",
" ... \\\n",
" 164834876 164834877 164834878 ... 164835261 164835262 164835263 \n",
"base T A T ... T T A \n",
"\n",
" \n",
" 164835264 164835265 164835266 164835267 164835268 164835269 164835270 \n",
"base T T A A T A A \n",
"\n",
"[1 rows x 353 columns]"
]
},
"execution_count": 83,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Obtain the molecule consensus sequence as pandas df:\n",
"pd.DataFrame({'base':molecule.get_consensus()}).T"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Visualisation"
]
},
{
"cell_type": "code",
"execution_count": 84,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<h3>chr1:164834869-164835270\n",
" sample:APKS3-P19-1-1_115 valid molecule</h3>\n",
" <h5>UMI:GAC Mapping qual:60.0 Cut loc: chr1:164834865 </h5>\n",
" <div style=\"white-space:nowrap; font-family:monospace; color:#888\"><h5>ST-E00285:221:HKM7VCCXY:6:1224:9120:30439</h5>AGTTAGATATGGACTCTTCTTCAGACACTTTGTTTAAATTTTAAATTTTTTTCTGATTGCAAATTAGTATTTAATGTAGAAAATTTGGAGAATTCAGGAAAACCACAAAGAAGAAAATTTACTGTTAACTTTGTGC.........................................................................................................................................................................................................................................................................<br />...........................................................................................................................................TAGCCTTCTAGCACTTTATGTATAAGTGTGTTTTTTATGTGGGTGTAAGAAAGTTACATTTCCTTTTACCTACCTATATTACTAAAAATGTGTTATGAATATTTTCCATATCATTAAACATTCTTCTCAAGCATAACTTTAAATAACTGCA...............................................................................................................<br /><h5>NS500413:404:HJ3KHBGX5:2:22312:23759:6537</h5>AGTTAGATATGGACTCTTCTTCAGACACTTTGTTTAAATTTTAAATTTTTTTCTGATTGC.....................................................................................................................................................................................................................................................................................................................................................<br />..........................................................................................................................................................................................................................................................................................................................................AATTTATA<text style=\"color:red; font-weight:500\" >T</text>TATAGAA<text style=\"color:red; font-weight:500\" >A</text>ATTT<text style=\"color:red; font-weight:500\" >T</text>A<text style=\"color:red; font-weight:500\" >T</text><text style=\"color:red; font-weight:500\" >T</text>TACTTTTG<text style=\"color:red; font-weight:500\" >T</text>TTTT..TTTTTTTTTTTTTTTTTTACTATTATTAATA<br /><h5>NS500413:404:HJ3KHBGX5:4:13507:22698:18344</h5>AGTTAGATATGGACTCTTCTTCAGACACTTTGTTTAAATTTTAAATTTTTTTCTGATTGC.....................................................................................................................................................................................................................................................................................................................................................<br />.................................................................................................................................................................................................................................................................................................................................................................................GTTTTTTTTTTTTTTTTTTACTATTATTAATA<br /><h5>NS500413:404:HJ3KHBGX5:3:23502:18585:5263</h5>AGTTAGATATGGACTCTTCTTCAGACACTTTGTTTAAATTTTAAATTTTTTTCTGATTGC.....................................................................................................................................................................................................................................................................................................................................................<br />...................................................................................................................................................................................................................................................................................................................................................................................TTTTTTTTTTTTTTTTTACTATTATTAATA<br />AGTTAGATATGGACTCTTCTTCAGACACTTTGTTTAAATTTTAAATTTTTTTCTGATTGCAAATTAGTATTTAATGTAGAAAATTTGGAGAATTCAGGAAAACCACAAAGAAGAAAATTTACTGTTAACTTTGTGC...TAGCCTTCTAGCACTTTATGTATAAGTGTGTTTTTTATGTGGGTGTAAGAAAGTTACATTTCCTTTTACCTACCTATATTACTAAAAATGTGTTATGAATATTTTCCATATCATTAAACATTCTTCTCAAGCATAACTTTAAATAA.............................................AATTTATA<text style=\"color:red; font-weight:800\" >T</text>TATAGAA<text style=\"color:red; font-weight:800\" >A</text>ATTT<text style=\"color:red; font-weight:800\" >T</text>A<text style=\"color:red; font-weight:800\" >T</text><text style=\"color:red; font-weight:800\" >T</text>TACTTTTG<text style=\"color:red; font-weight:800\" >T</text>TTTT.GTTTTTTTTTTTTTTTTTTACTATTATTAATA<br />AGTTAGATATGGACTCTTCTTCAGACACTTTGTTTAAATTTTAAATTTTTTTCTGATTGCAAATTAGTATTTAATGTAGAAAATTTGGAGAATTCAGGAAAACCACAAAGAAGAAAATTTACTGTTAACTTTGTGC???TAGCCTTCTAGCACTTTATGTATAAGTGTGTTTTTTATGTGGGTGTAAGAAAGTTACATTTCCTTTTACCTACCTATATTACTAAAAATGTGTTATGAATATTTTCCATATCATTAAACATTCTTCTCAAGCATAACTTTAAATAA?????????????????????????????????????????????AATTTATA<text style=\"color:black; font-weight:800\" >C</text>TATAGAA<text style=\"color:black; font-weight:800\" >C</text>ATTT<text style=\"color:black; font-weight:800\" >A</text>A<text style=\"color:black; font-weight:800\" >G</text><text style=\"color:black; font-weight:800\" >C</text>TACTTTTG<text style=\"color:black; font-weight:800\" >A</text>TTTT?GTTTTTTTTTTTTTTTTTTACTATTATTAATA<br /></div>"
],
"text/plain": [
"<IPython.core.display.HTML object>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"# Display the molecule here in the notebook:\n",
"from IPython.core.display import display, HTML\n",
"display(HTML( molecule.get_html() ))"
]
},
{
"cell_type": "code",
"execution_count": 85,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"AGTTAGATATGGACTCTTCTTCAGACACTTTGTTTAAATTTTAAATTTTTTTCTGATTGCAAATTAGTATTTAATGTAGAAAATTTGGAGAATTCAGGAAAACCACAAAGAAGAAAATTTACTGTTAACTTTGTGC........................................................................................................................................................................................................................................................................."
],
"text/plain": [
"<IPython.core.display.HTML object>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"# Display a single read:\n",
"fragment = molecule[0]\n",
"display(HTML(fragment.get_html(span_start=molecule.spanStart, span_end=molecule.spanEnd,show_read1=1,show_read2=0) ))"
]
},
{
"cell_type": "code",
"execution_count": 86,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"\"ST-E00285:221:HKM7VCCXY:6:1224:9120:30439\\t99\\t0\\t164834865\\t60\\t140M\\t0\\t164835008\\t140\\tCATGAGTTAGATATGGACTCTTCTTCAGACACTTTGTTTAAATTTTAAATTTTTTTCTGATTGCAAATTAGTATTTAATGTAGAAAATTTGGAGAATTCAGGAAAACCACAAAGAAGAAAATTTACTGTTAACTTTGTGC\\tarray('B', [32, 37, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 37, 41, 41, 41, 41, 41, 37, 37, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 41, 37, 41, 37, 32])\\t[('NM', 0), ('MD', '140'), ('MC', '151M'), ('AS', 140), ('XS', 20), ('Is', 'ST-E00285'), ('RN', '221'), ('Fc', 'HKM7VCCXY'), ('La', '6'), ('Ti', '1224'), ('CX', '9120'), ('CY', '30439'), ('Fi', 'N'), ('CN', '0'), ('aA', 'TGACCAAT'), ('LY', 'APKS3-P19-1-1'), ('RX', 'GAC'), ('RQ', 'AAF'), ('BI', '115'), ('bc', 'GCACACGC'), ('BC', 'GCACACGC'), ('QT', 'FFJJJJJJ'), ('MX', 'NLAIII384C8U3'), ('MI', 'GCACACGCGACTGACCAAT'), ('QM', 'FFJJJJJJAAF////////'), ('SM', 'APKS3-P19-1-1_115'), ('rS', 'TGCAGT'), ('rP', 164835159), ('RC', 1), ('DT', 'NLA'), ('RG', 'HKM7VCCXY.6.APKS3-P19-1-1_115'), ('DS', 164834865), ('RS', 1), ('RZ', 'CATG')]\""
]
},
"execution_count": 86,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"str(fragment[0])"
]
},
{
"cell_type": "code",
"execution_count": 88,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]\n",
"[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]\n",
"[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]\n",
"[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]\n",
"[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]\n",
"[4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]\n",
"[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]\n",
"[4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]\n"
]
},
{
"data": {
"text/plain": [
"<matplotlib.image.AxesImage at 0x7f0c29101c50>"
]
},
"execution_count": 88,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAfMAAANjCAYAAAC+9vMhAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAmcgAAJnIBISuX7gAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xm0bWdZJ+rfS3JIa2gEDDcJjUUTUBIa0QESmxooIFIiioiUpcYGqVKxoqXWEC0FpUSLeMvSogApqLIUO64CcqWtK7FBGkWwBCIigSQGEgwQEg5p3/vHWsczz+bsc9Y+e+291rfP84yxxpjzW9/85rt3zspvz7nm/GZ1dwCAcd1u1QUAANsjzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAa3Z8O8qh5ZVS+qqvdU1XXz13vmbY9cdX0AsCzV3auuYamq6rQkv5TkwqN0/e9JfqC7b9j5qgBg5+ypMK+qE5L8v0m+etK8P8nfJLklyQOTnDF573VJHt/dt+5akUmq6oYkJyW5LcnVu7lvAHbc3TI7831jd5+2Gzvca2H+3CT/ftL04iQ/1t3Xzt8/LcmPJXnWpM9zu/vHd6/KpKpuSXLCbu4TgF13a3efuBs72jNhXlVnJfm7JCfPm36tu//VJn2fk4OBvj/Jfbr7H3a+yn/a/01J9iXJSTllt3bLsu1b4md0mR/DWuJYyxxsmf+vqWX+kGv6/8BllbXM39VS/xsub6h1+/zcePP1mRd1c3fffvsjHt1eCvPnJfmR+eqnk5xz4Ij8MH1vn1nwnzNv+vnu/tGdr/Kf9n9FkrNOyim5oB6/W7tlyW533rlLG6tuuW1pY912yr6ljdUnLO8a2ROuv3FpY916+klLG6tuXd7vfplut//mpYyzzN/VMv8bLvPf6bJ+V8ly6rrkXRfnxps/lSRXdvfZ2x5wAXvpavYnTZZ/e7MgT5LuvinJSydNX79jVQHADtsTYV5V909yn0nTaxfY7A8ny/etqvsttyoA2B17IsyTnL9h/S0LbPOXSW46whgAMIS9EuYPmCzflOTyo20wP9U+7feAzfoCwDrbK2F+z8nyFb34VX0fnizfa3nlAMDu2ZX733bBdCKYT25hu+smy5+znQKq6qIkFy3Y/czt7AsApvZKmE9n2PnMFrbbv8kYx+KMJGdtcwwA2LK9EubTGwNv2cJ2077bvbH/uiRXLtj3zJgBDoAl2Sth/unJ8smb9vps077beuBKd1+c5OJF+h6YNGY7+wOAA/bKBXDXT5a3Mj/qqZuMAQDD2Cth/rHJ8t23sN30QrR/XFItALCr9kqYXzpZ/tyqOnXTnoc6Z7L8viXWAwC7Zq+E+Xs3rD/4aBvMn7J21yOMAQBD2Cth/rYk08f5PGqBbS6YLH9mPgYADGdPhHl3X5/kTZOmpy2w2bTPm7p7W1ezA8Cq7Ikwn3vZZPm8qnrCZh2r6qFJHrfJtgAwlL1yn3mS/G6Sd+Xg089eWFXv7+5DLmyrqrsn+V85OGnLXyV5xa5VyZ5x27tdM7kVt666gE0s+iCH3baOv691rClZw7p6K3OXLceeCfPu7qr6riSXZHav+d2TvLWqXjBvuyXJFyf5viSfN99sf5Lv3sKDWQBg7eyZME+S7n5HVT0tya9nFuhnJPnR+Wuj/Ume1t3v2MUSAWDp9tJ35kmS7v69JA/L7IK4wx1xd5I3JnnovC8ADG1PHZkf0N3vTfLoqjonySNzcB70K5P8WXdfvrLiAGDJ9mSYHzAP7d9adR0AsJP23Gl2ADjeCHMAGJwwB4DBCXMAGJwwB4DBCXMAGJwwB4DBCXMAGJwwB4DBCXMAGJwwB4DBCXMAGJwwB4DB7emnpsFOut155y5trLrltqWNddsp+5Y2Vp+wvL/3T7j+xqWNdevpJy1trLp1eb/7Zbrd/puXMs4yf1fL/G+4zH+ny/pdJUuq612vS5ZX0kIcmQPA4IQ5AAxOmAPA4IQ5AAxOmAPA4IQ5AAxOmAPA4IQ5AAxOmAPA4IQ5AAxOmAPA4IQ5AAxOmAPA4IQ5AAxOmAPA4IQ5AAxOmAPA4IQ5AAzuxFUXAKO67d3vW3UJQ7l11QVsolddwCbW8fe1jjUla1hX37Lru3RkDgCDE+YAMDhhDgCDE+YAMDhhDgCDE+YAMDhhDgCDE+YAMDhhDgCDE+YAMDhhDgCDE+YAMDhhDgCDE+YAMDhhDgCDE+YAMDhhDgCDO3HVBcCobnfeuUsbq265bWlj3XbKvqWN1Scs7+/9E66/cWlj3Xr6SUsbq25d3u9+mW63/+aljLPM39Uy/xsu89/psn5XyZLqetfrkuWVtBBH5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIPbU2FeVV9RVX0Mr3NXXTsAHKs9FeYAcDw6cdUF7KDPJHnzgn2v38lC2Jtue/f7Vl3CUG5ddQGb6FUXsIl1/H2tY03JGtbVt+z6LvdymH+0ux+76iIAYKc5zQ4AgxPmADA4YQ4AgxPmADA4YQ4Ag9vLYX7HqvrtqrqsqvZX1aeq6oNV9ftV9X1VdcaqCwSAZdjLYX6HJE9Ocs8kJyc5Pcm9knxdkv+S5MNV9f0rqw4AlmQv32eeJJcluTLJjUnukuSBOfgz3yHJL1XVg7v7O7e7o6q6KMlFC3Y/c7v7A4AD9lqY35bkjUlemuS13X3t9M2qOj3Jv0zynMzCPUkurKq/7e7nbXPfZyQ5a5tjAMCW7akw7+5LknzVEd6/Psl/q6rXJLkks9PuSfKTVfWy7v7oNnZ/XWZnARZxZpITtrEvAPgne/k780119+VJvnnSdGqSbZ1q7+6Lu/vsRV5JPrKdfQHA1HEZ5knS3W9N8keTpk2P6AFgnR23YT73R5Pl+62qCADYjuM9zK+aLN9l014AsMaO9zA/dbL86ZVVAQDbcLyH+QMny1evrAoA2IbjNsyr6pQk/2LS9GerqgUAtmNP3We+Rc9J8nmT9d9fVSGM6Xbnnbu0seqW25Y21m2n7FvaWH3C8v7eP+H6G5c21q2nn7S0serW5f3ul+l2+29eyjjL/F0t87/hMv+dLut3lSyprne9LlleSQvZM0fmVfXVVfX8qjr7KP32VdXPJfmhSfNfJnnVjhYIADtkLx2Zn5rZ3Og/WFV/muTNSf5Pko8luSmzq9W/OMnTkpwz2e7aJN/S3b275QLAcuylMD/gdkkumL+O5v1JntLdl+5sSQCwc/bMafYk78vse++PL9D3siQ/kuQh3f3OnSwKAHbanjky7+73Jfn6JKmqf5bkAUnOTnLHzB5qcl1mt5+9vbv/flV1AsCy7Zkwn+ruDyT5wKrrAIDdsJdOswPAcUmYA8DghDkADE6YA8DghDkADE6YA8DghDkADE6YA8DghDkADE6YA8DghDkADE6YA8Dg9uSDVmA33Pbu9626hKHcuuoCNtGrLmAT6/j7WseakjWsq2/Z9V06MgeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwa00zKvqrlX1uKr6yap6VVVdVVU9eX37Nsb+/Kp6dlX9RVVdU1X7q+oDVfV7VfWNVXXCEn8UAFiZE1ex06o6M8mfJ7nnDo3/g0l+LslJG976/PnriUneUlVP6+4P7kQNALBbVnVkfnJ2Lsh/Iskv5mCQ35bk/yS5JMlVk66PSHLJ/A8LABjWOnxnfk2S1yb5mSRft52BquoxSX560vSWJA/o7gd195cnOTvJU5NcP3//7CS/u519AsCqreQ0e5Jrkzw5ydu7+0PTN6rqmAas2YbPS3JggEuTPLq7P32gT3ffluQ3q+pjSd4wb/7Sqvr67v69Y9oxAKzYSo7Mu/u67v7djUG+TY9Ncv5k/ZnTIN+w/zcm+a1J048tsQ4A2FXrcJp9Wb5hsvzBJK8/Sv8XTpYfXlVnL78kANh5eynMHz9Zfl1391H6/3GSG+bLleRrdqQqANhheyLMq+puSaZXpb/laNt09y1J3j5pOn+zvgCwzvZEmCd5wIb1Dyy43bTfxjEAYAh7Jcw33rP+4QW3m/a713JKAYDdtapb05btjA3rn1xwu+smy5+znQKq6qIkFy3Y3UQ1ACzNXgnz0zasf2bB7fYfYYytOiPJWdscAwC2bK+E+b4N67csuN203+23WcN1Sa5csO+ZSTzoBYCl2CthvnFymJMP03Y4J0+Wb9i01wK6++IkFy/St6quiKN4AJZkr1wAd/2G9VMW3O7UI4wBAEPYK2H+sQ3rd19wu+mFaP+4pFoAYFftlTC/dMP6PRbc7pzJ8vuWVAsA7Kq9Eubvz6EXsz14we0eMll+7/LKAYDdsyfCvLtvSvLWSdOjjrZNVZ2Z5D6TpkuWXRcA7IY9EeZzr5wsP3o+X/uRPG2y/Ikkb15+SQCw8/ZSmL88yY3z5X1JfnSzjlV1epIfmDT9+vzBKwAwnD0T5t19RQ59Rvkzq+obNvarqn1JXpqDF8ntT/KzO18hAOyMlYV5Vb24qj6z8bWh2yJ9pn4qs4vhktkMa79dVb9WVd9QVV9ZVd+b5B1JvnGyzb/r7quW9oMBwC5b5Qxw+5KcdJQ+J2YLNXb3x6vqa5O8MbPbzm6X5F/OX4fzvO7+lUXHB4B1tGdOsx/Q3X+b5LwkL8mhD1KZek+SJ3T3j+1aYQCwQ1Z2ZN7d357k23do7E8k+a6q+rdJ/nlmR+mnJbkqyV939zt3Yr8AsAp75UErh9Xdn8qht6wBwJ6z506zA8DxRpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMbqVhXlV3rarHVdVPVtWrquqqqurJ69u3OF4fw+t7d+jHA4BdceIqdlpVZyb58yT3XMX+AWAvWUmYJzk5Ox/klyTZv0C/D+9wHQCwo1YV5lPXJPmLJO9I8vYkr1zSuN/W3ZctaSwAWFurCvNrkzw5ydu7+0PTN6pqNRUBwKBWEubdfV2S313FvgFgr3FrGgAMTpgDwOCEOQAMbh2uZt8pv1BVD0xyTpJ9Sf4xyfuTvDnJ/+juD66yOABYlr0c5t+4Yf2s+esrkvx4Vb0kyb/t7kXuRT+qqrooyUULdj9zGfsEgGRvh/k1Sf4+yfVJ7pDk3CSnz987McnTk3xxVX1ld39yCfs7I7M/FgBgV+21MH9PkhcleXV3//30jao6Mcljkjw3yXnz5ock+c0kj1vCvq9LcuWCfc9McsIS9gkAeyvMu/sLjvDeLUleU1Vvyuwe98fP33psVT2hu1+9zX1fnOTiRfpW1RVxFA/Akhx3V7N392eSPDXJRyfN37+icgBg2467ME+S7v5UkhdMmi6oqpNXVQ8AbMdxGeZzfzRZPjmzW9gAYDjHc5hftWH9LiupAgC26XgO81M3rH96JVUAwDYdz2H+wA3rV6+kCgDYpuM5zL95snxZd2887Q4AQzguw7yqnpDkCZOm319VLQCwXXsizKvqDlX1iqp62AJ9n5rk5ZOmTyd53o4VBwA7bGUzwFXVi5N861G6vbiq/tvGxu7eeE94JXlSkidV1fuSvC7JX2V2xfoNST4nyYMye/jKw6dDJfmO7v7IMf0QALAGVjmd674kJx2lz4nZeo3nzl9H86kkT+/u397i+ACwVvbEafYk+5O8NMmHFuj7ySS/lOQLu/vlR+sMAOtuZUfm3f3tSb59SWPdmOTCJKmqu2X2VLR7ZDYRzEmZPQb12iTvTvLu7r51GfsFgHWwp56aliTdfXWSN666DgDYLXvlNDsAHLeEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAMTpgDwOCEOQAM7sRV7biq7pjkK+evBye5X5I7Jbk5ybVJ3pXkTUn+R3d//BjGf1CS70jy6CRnJ7l9kiuTvCPJr3X3a5fwYwDAyu16mFfVuUl+IclXZxawG90+yWlJzknytUl+tqp+PMl/7u5eYPwTkzw7yY/ms8883G/++paq+oMkF3b3Ncf6swDAOljFafYvzCykp0F+a5JLk1yS5E8zOzI/4NQkv5jkxVVVC4z/wiT/Pgd/tpszO8r/0yT/OOn3tUneUFWnHcPPAABrY5Xfmd+S5PeTPDHJnbv73O7+8u5+VJK7zNuvnPT/ziTfe6QBq+p7klw4aXpVknt394Pn4949yffP950k5yd50TJ+GABYlVWE+c1JfjXJP+vur+/uV3b3ddMOPfPKJI9I8pHJW8+uqn2HG7SqTk3y05OmP0rypO7+pz8Iuvvm7v7lHPpHwVOr6qHb+okAYIV2Pczn4f3d3f3hBfpenuQ/TJrukuTLNun+bUnOPLBpkmd0962bjPuSJG+dr1Zm368DwJBGuDXt1RvWz92k3zdMlt/c3e87yrgvnCx/TVWdtOXKAGANjBDm125YP2Njh6o6PYcesS9y29kfTpZPT/LlWy8NAFZvhDC/54b1qw/T54FJpt+lv+Vog3b3R5JcNmk6f8uVAcAaGCHMn7Rh/c8P0+cBG9Y/sODY034bxwCAIax1mFfVHZI8c9L07u7+m8N0nR6935LkqgV3Mb0I715bqw4A1sPKpnNd0PNz8Ar1JHnWJv2m36N/qrtvW3D86S1xn7OVwjaqqouSXLRg9zOP3gUAFrO2YV5VF2Y2UcwBv9XdG69sP2A6i9tntrCb/ZuMcSzOSHLWNscAgC1byzCvqguS/NdJ0weTPP0Im0wvfrtl016fbdr3cPPEb8V1OXTGuiM5M8kJ29wfACRZwzCvqvMzu7f8wH3fVyd5bHd/8gibfXqyfPIWdjfte8MWtvss3X1xkosX6VtVV8RRPABLslYXwFXV/ZO8Pskd5k0fT/LV3f23R9n0+snyKVvY5ambjAEAw1ibMK+qeyd5Y5K7zZuuT/K47n7XApt/bLJ8+nwSmUVML0T7x017AcAaW4swr6qzk7wpydnzpv1Jvra737r5Voe4dMP6PRbc7pzJ8tGmfwWAtbTyMK+qz8vsiPze86Ybkzyxu9+8hWHeu2H9wQvsd19mz1bfbAwAGMJKw7yq7pzkDUnuP2+6Ock3dffrtzJOd/99kismTY9aYLOH5dDv1y/Zyj4BYF2sLMyr6owkr0vyoHnTrUme1t2vOsYhp9s9uaqOdqvZ0ybLf9Pdi04BCwBrZSVhXlWnJnlNki+aN92W5Du6+3e2MezLJst3yRHuS59/R/9tm2wLAEPZ9TCfPzf8lTl4KryTfE93/9p2xu3ut+fQo/PnVtWXHmb/ZyT5jRycvvWqJL+ynX0DwCqtYtKYZyZ59GT9E5mdFn/ygtu/obuff4SxH5Hkrpk9o/xNVfWSzO5dvz7JeUm+PwcvtrstydO7e/9hxgKAIawizE/dsH6nJI/ZwvYf2eyN7r6sqp6Y2Qxyd85sFrl/PX9tdGuSZx5hvncAGMLKb01btu7+s8yOwF+Rzedpf1uSC7rb6XUAhrfrR+bd/VNJfmqH93Flkm+sqrsm+bLMJqO5fZJ/SPL2BaaHBYBhrN2DVpapu6/J7AgdAPasPXeaHQCON8IcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAYnzAFgcMIcAAZ34qoLgFE94/1/t7Sxzr/9R5Y21r33nb60sYCtu8dDb8qVV+3uPh2ZA8DghDkADE6YA8DghDkADE6YA8DghDkADE6YA8DghDkADE6YA8DghDkADE6YA8DghDkADE6YA8DghDkADE6YA8DgVvY886q6Y5KvnL8enOR+Se6U5OYk1yZ5V5I3Jfkf3f3xBca7V5IPHkMpj+vu1x7DdgCwFnY9zKvq3CS/kOSrk9z+MF1un+S0JOck+dokP1tVP57kP3d371qhADCIVRyZf2FmIT11a5K/S/LRJCckeUCSO8/fOzXJLyb5wqr67i0E+usW7HfNgv0AYC2t7DR7kluS/EGSlyX5/7r7ugNvVFUl+RdJfiXJWfPm70zyF0lesMjg3f3YZRYLG73gvvdZ4mjLHAtYpY/1lUn27+o+V3EB3M1JfjXJP+vur+/uV06DPEl65pVJHpHkI5O3nl1V+3axVgBYe7se5vPw/u7u/vACfS9P8h8mTXdJ8mU7VhwADGiEW9NevWH93JVUAQBraoQwv3bD+hkrqQIA1tQIYX7PDetXr6QKAFhTI4T5kzas//kiG1XV/6yq91fVDfPXh6vqtVX1I1V1tx2oEwBWYq3DvKrukOSZk6Z3d/ffLLj5t2Z2v8+p89c5SR6T5HlJPlRVz6mqE5ZZLwCswirvM1/E85OcOVl/1ha2vSrJZZnd7HenzCaiOXn+3snzsR5eVU/o7pu3W2hVXZTkogW7n3n0LgCwmLUN86q6MLOJYg74re7eeGX7VCd5W5IXJ3lNd1+1YbyTMjtl/9wk95o3PybJLyV5xhJKPiMHJ7gBgF2zlmFeVRck+a+Tpg8mefqRtunuDyX5kiO8f2OSl1fVa5O8IcnD5m89vape0N3v3l7VuS7JlQv2PTOzaWsBYNvWLsyr6vzM7i0/ad50dZLHdvcnlzF+d3+8qp6U5NLMTrdXku9L8j3bHPfiJBcv0reqroijeACWZK0ugKuq+yd5fZI7zJs+nuSru/tvl7mf+exzvzlp+qpljg8Au2ltwryq7p3kjUkO3DZ2fWbPGn/XDu3yjybL96qqwz2OFQDW3lqEeVWdneRNSc6eN+1P8rXd/dYd3O1VG9Y/dwf3BQA7ZuVhXlWfl9kR+b3nTTcmeWJ3v3mHd33qhvVP7/D+AGBHrDTMq+rOmV1Zfv95081Jvqm7X78Lu3/gZPnGZV1gBwC7bWVhXlVnJHldkgfNm25N8rTuftUu7LuSPGXS9Jad3icA7JSVhHlVnZrkNUm+aN50W5Lv6O7f2aUS/k2S8ybrv79L+wWApdv1+8znM7G9Msmj5k2d5Hu6+9e2MeYXZDaV6s9396VH6FdJfiCzaWIP+IckLzrWfXP8esb7/25pY9133zVLGwtYrUd/yY25+iO7u89VTBrzzCSPnqx/IsmTq+rJC27/hu5+/oa2fUkuTHJhVf1Fkv+d5F2ZTThzYG72hyR5apJzJ9vdmOQp3b1/yz8FAKyJVYT5xqvI75TZHOmLOtrfOw/LwalajzbOt3b3n2xh3wCwdlZ+a9qSXJXk5Tl60CfJR5P8TJIHdfcbd7QqANgFu35k3t0/leSnljzmR5N8S5JU1TlJviDJPTI76t+X5FNJPpbknUne2929zP0DwCqt3YNWtqu7L09y+arrAIDdsldOswPAcUuYA8DghDkADE6YA8DghDkADE6YA8DghDkADE6YA8DghDkADE6YA8DghDkADE6YA8DghDkADG7PPTUNdssL7nufJY62zLGAVfpkvybJ/l3dpyNzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABicMAeAwQlzABjciasuAEb1jPf/3apLOKz77rtm1SUM5Zvf+Z1LG+s3H/KSpY3FuB79JTfm6o/s7j4dmQPA4IQ5AAxOmAPA4IQ5AAxOmAPA4IQ5AAxOmAPA4IQ5AAxOmAPA4IQ5AAxOmAPA4IQ5AAxOmAPA4IQ5AAxOmAPA4IQ5AAxOmAPA4IQ5AAyuunvVNRx3quqKJGedlFNyQT1+1eUAsER/3K/JjdmfJFd299m7sU9H5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAINFXUBLAAATP0lEQVQT5gAwOGEOAIM7cRU7rap9SR6a5JFJzktybpJ7Jrljkn1JrktyeZJ3JPmdJG/o7tu2uI8HJfmOJI9OcnaS2ye5cj7mr3X3a5fywwDAiq0kzJM8N8kPH+H9O89f5yf5ziR/VVUXdvc7jzZwVZ2Y5NlJfjSffebhfvPXt1TVHyS5sLuvOYb6AWBtrCrMa8P6DUk+kOTjSTrJmZmF7oEwfnCSS6rqsd39p0cZ+4VJLpys35zkPUmuz+wMwOfO2782yRuq6ku7+4Zj/UE4fj3j/X+3tLHuu8/flKvy/pvvurSx/Hdcjf9w+ROWNtZPn/PqbY/x6C+5MVd/ZAnFbMGqvjPfn+QPknxPknO7+/TuPr+7v6K7v7K7H5BZoP9sklvn25ye5Deq6rTNBq2q78mhQf6qJPfu7gd396OS3D3J9ye5Zf7++UletMwfDAB220qOzLv7Jxboc02SZ1XVZUlePG++R5JvSvLSjf2r6tQkPz1p+qMkT+ruA38MpLtvTvLLVbU/ya/Om59aVc/v7r88hh8FAFZu7a9m7+5fzewU/AFfsUnXb8vsaD6Znap/xjTIN4z5kiRvna9WZt+vA8CQ1j7M56ZHzWdu0ucbJstv7u73HWXMF06Wv6aqTjqmygBgxUYJ8+nXAddtfLOqTk/yZZOmRW47+8PJ8ulJvvzYSgOA1Vr7MJ/fk/6ISdNbDtPtgZndn36kPofo7o8kuWzSdP6x1AcAq7b2YZ7ZFe0HTq1fm+Rlh+nzgA3rHzhMn8OZ9ts4BgAMYVX3mW9qPunLXZN8SZJ/neSr5m99Jsm3dPe1h9nsnpPlW5JcteDuPjxZvtfWKj1UVV2U5KIFu2/2vT8AbNlahHlVfSwHJ3M5nDcm+aHufvcm758xWf7UFqZ+nX7//jkLbrOZM5Kctc0xAGDL1iLMj+JPkvxykr8+Qp/pRDKf2cLY+zcZ41hcl9nc74s4M8kJ29wfACRZnzB/U5I7zJdPyqHTuT5q/npbVT2luy87zPbTi99uOcz7m5n2vf0Wtvss3X1xkosX6VtVV8RRPABLshZh3t1P2dhWVXdO8l1JfjKzo+YvTvLmqnp4d1+9ofunJ8snb2HX077mZwdgSGt7NXt3X9vdP5/kgiSfmjffI8nzD9P9+snyKVvYzambjAEAw1jbMD9g/tjT/zhp+ub5UfvUxybLp88nkVnE9KryfzyW+gBg1dY+zOd+e7J8YpIv2vD+pRvW77HguOdMlo82/SsArKVRwvzyDet32bD+3g3rDz7agPOZ5b7wCGMAwBBGCfM7bFj/xHSlu/8+yRWTpkctMObDcuj365ccW2kAsFqjhPkFG9YPN13rqybLT66qo91q9rTJ8t9096JTwALAWqnuXnUNRzQP5T9P8pB50we6+z6H6ffwJG+bNP1Ad/+XTcY8O8l7cnDWt3/X3f9peVUf2YH7zE/KKbmgHr9buwVgF/xxvyY3zuYku7K7z96Nfe76kXlVfVVV/UJVHXXSlKq6e5JX52CQJ8nPHa5vd789hx6dP7eqvvQwY56R5DdyMMivSvIrC5YPAGtnFZPGnJbkh5P8UFX9WZI/zmyq1msym/zl9CSfn9mp9a/LofeCvyrJS44w9jMze1zqXefjvKmqXpLk9ZndR35eku9Pcu95/9uSPL279x9mLAAYwipngKskXzp/LeKlSb63j/C9QHdfVlVPzOxo/s6ZTQ37r+evjW5N8szufvWWqgaANbOKC+Dekdkc5u9JcrQv7G9K8ookX97dF3b3TUcbvLv/LLMj8Fdk83na35bkgu52eh2A4e36kXl3X5HkhzI7zX7HJOdndlr9LpkdSd+Q5OOZ3ff9ru7eylPQDuzjyiTfWFV3TfJlSc7O7EEq/5Dk7d39t8v4WQBgHaz0QSvd/Ykkb56/dmL8azI7QgeAPWuU+8wBgE0IcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAYnDAHgMEJcwAY3Imr2GlV7Uvy0CSPTHJeknOT3DPJHZPsS3JdksuTvCPJ7yR5Q3ffdpQx75Xkg8dQzuO6+7XHsB0ArIWVhHmS5yb54SO8f+f56/wk35nkr6rqwu5+524UBwAjWVWY14b1G5J8IMnHk3SSM5PcLwe/Bnhwkkuq6rHd/acL7uN1C/a7ZsF+ALCWVhXm+5P8QZJXJbmkuy/d2KGq7prkmUl+LMkJSU5P8htV9cDuvuFoO+juxy63ZABYTysJ8+7+iQX6XJPkWVV1WZIXz5vvkeSbkrx056oDgLGs/dXs3f2rmZ2CP+ArVlQKAKyltQ/zub+cLJ+5sioAYA2NEubTrwOuW1kVALCG1j7M5/ekP2LS9JZV1QIA62jtwzzJz+bgqfVrk7xskY2q6n9W1fur6ob568NV9dqq+pGquttOFQsAu23twryqTqyqu1fVE6vq9Un+3fytzyT5lu6+dsGhvjXJfZKcOn+dk+QxSZ6X5ENV9ZyqOmHJ5QPArlvVfeaHqKqPJfncI3R5Y5If6u53b2HYq5Jcltk97XdK8oAkJ8/fOznJs5I8vKqe0N03b7noDarqoiQXLdjdRXwALM1ahPlR/EmSX07y10fp10neltk96a/p7qumb1bVSUmelNlUsveaNz8myS8lecYS6jwjyVlLGAcAtmRdwvxNSe4wXz4ph07n+qj5621V9ZTuvuxwA3T3h5J8yWY76O4bk7y8ql6b5A1JHjZ/6+lV9YItHvUfznVJrlyw75mZzWoHANtW3b3qGg6rqu6c5LuS/GSS0+bNH07y8O6+eptj3yPJpTl42v3F3f092xlzi/u/IslZJ+WUXFCP363dArAL/rhfkxuzP0mu7O6zd2Ofa3cB3AHdfW13/3ySC5J8at58jyTPX8LYH07ym5Omr9rumACwKmsb5gfMH3v6HydN3zw/at+uP5os36uqbr+EMQFg1619mM/99mT5xCRftIQxr9qwfqSr6QFgbY0S5pdvWL/LEsY8dcP6p5cwJgDsulHC/A4b1j+xhDEfOFm+sbs/uYQxAWDXjRLmF2xY/8Bhey2oqirJUyZN5nsHYFhrH+bzC9OeNWn6QHdfus1h/02S8ybrv7/N8QBgZXY9zKvqq6rqF6rqqLOlVdXdk7w6yUMmzT93mH5fUFUvqar7H2W8qqpnJvm/J83/kORFi1UPAOtnFTPAnZbkh5P8UFX9WZI/zmyq1msyuwjt9CSfn9mp9a/LoReqvSrJSw4z5r4kFya5sKr+Isn/TvKuJFfn4NzsD0ny1CTnTra7MclTunv/sn44ANhtq5zOtZJ86fy1iJcm+d4++pR1D8vBqVqP5CNJvrW7/2TB/QPAWlrFd+bvSHJxkvdk9nCUI7kpySuSfHl3X9jdN23S76okL88soI/mo0l+JsmDuvuNi5UMAOtr14/Mu/uKJD+U2Wn2OyY5P7PT6nfJ7CErNyT5eJL3JnlXd39mgTE/muRbkqSqzknyBZlN/XqnzE7BfyrJx5K8M8l7Fzi6B4BhrPSpad39iSRvnr+WNebl+exJZgBgz1r7W9MAgCMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwuOruVddw3Kmqm5LsS5KTcsqKqwFgmW7M/gOLN3f37Xdjnyfuxk74LP90RmTyHx2AvWXXzn4L89W4MclJSW5LcvUmfc5MckKSW5N8ZJfqglXyb5694m6ZBfmNu7VDp9nXVFVdkeSsJFd299mrrgd2mn/zcOxcAAcAgxPmADA4YQ4AgxPmADA4YQ4AgxPmADA4YQ4AgxPmADA4YQ4AgzOd6/q6OMkZSa5bdSGwS/ybh2NkOlcAGJzT7AAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5mukqh5ZVS+qqvdU1XXz13vmbY9cdX1wJFV116p6XFX9ZFW9qqquqqqevL59G2N/flU9u6r+oqquqar9VfWBqvq9qvrGqjphiT8KDMdT09ZAVZ2W5JeSXHiUrv89yQ909w07XxUspqrOTPLnSe55lK7f0d0vO4bxfzDJzyU56Qjd3pLkad39wa2OD3uBI/MVmx9R/D85NMj3J3lHZv+DnD7b+cIkr3AUwpo5OUcP8mNSVT+R5BdzMMhvS/J/klyS5KpJ10ckuWT+hwUcd4T56j0nyVdP1l+c5Ozufnh3PyLJ/5XkZybvPybJs3exPtiKa5K8NrN/s1+3nYGq6jFJfnrS9JYkD+juB3X3lyc5O8lTk1w/f//sJL+7nX3CqJxmX6GqOivJ32V2ZJMkv9bd/2qTvs9J8qz56v4k9+nuf9j5KuHIquqMzP4gfXt3f2jDe9P/wSx8mr2qKsk7k5w/b7o0yUO7+9OH6fvoJG+YND2pu39v8Z8AxufIfLV+IAeD/NNJfvAIfZ+T5PL58ilJnrmDdcHCuvu67v7djUG+TY/NwSBPkmceLsjn+39jkt+aNP3YEuuAIQjz1XrSZPm3u/vazTp2901JXjpp+vodqwpW7xsmyx9M8vqj9H/hZPnhVXX28kuC9SXMV6Sq7p/kPpOm1y6w2R9Olu9bVfdbblWwNh4/WX5dH/37wD9OcuAuj0ryNTtSFawpYb46529Yf8sC2/xlkpuOMAYMr6rulmR6VfpRPxvdfUuSt0+afDY4rgjz1XnAZPmmHPw+fFPzU+3Tfg/YrC8MbOO/6w8suN20n88GxxVhvjrT+3KvWOA04gEfnizfa3nlwNrYeM/6hw/b67P5bHDcEuarc8Zk+ZNb2G46icznLKkWWCdnbFhf9PPhs8FxS5ivzmmT5c9sYbv9m4wBe8XGf9eLfj58NjhuCfPV2TdZvmUL20373n5JtcA62bdhfdHPh88Gxy1hvjrTCTBO3rTXZ5v29cAV9qKNk8Ms+vnw2eC4JcxX5/rJ8ilb2O7UTcaAvWLjv+tFPx8+Gxy3hPnqfGyyfPctbDe9//Yfl1QLrJOPbVhf9PPhs8FxS5ivzqWT5c+tqlM37XmocybL71tiPbAuLt2wfo8Ft/PZ4LglzFfnvRvWH3y0DeZPWbvrEcaAveD9OfRitqN+NuYeMln22eC4IsxX521JbpysP2qBbS6YLH9mPgbsKfOZDt86aTrqZ6Oqzsyhzzq4ZNl1wToT5ivS3dcnedOk6WkLbDbt86budsUue9UrJ8uPns/XfiTTz8Ynkrx5+SXB+hLmq/WyyfJ5VfWEzTpW1UOTPG6TbWGveXkOnrnal+RHN+tYVacn+YFJ06/PH7wCxw1hvlq/m+Rdk/UXVtW5GztV1d2T/K8kJ8yb/irJK3a+PFiN7r4ihz6j/JlV9Q0b+1XVviQvzcGL5PYn+dmdrxDWSy3+fA92QlV9UWbf7x24l/a6JC+Yt92S5IuTfF+Sz5u/vz/Jl3X3O3a5VNhUVb04ybce5q2TJsu3JLl1Y4fuPuykMFV1p8y+O7/vvOm2JL+R5PeTXJvk/kmekeS8yWbf192/stX6YXTCfA1U1dcn+fUcfXKM/Ume1t2/t/NVweKq6mVJvu1Ytu3uOsK490vyxhx629lmntfdP3YsNcDonGZfA/NwflhmF8Qd7q+rzux/aA8V5BxPuvtvMzvyfkkOfZDK1HuSPEGQczxzZL5mquqcJI9Mcta86cokf9bdl6+uKli9qvqcJP88s6P005JcleSvu/udKy0M1oAwB4DBOc0OAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIMT5gAwOGEOAIP7/wFWeSmLD+gzYQAAAABJRU5ErkJggg==\n",
"text/plain": [
"<Figure size 1500x1000 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"tensor = molecule.get_alignment_tensor(\n",
" max_reads=8,\n",
" centroid=molecule.spanStart,\n",
" window_radius=10)\n",
"plt.imshow(tensor)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.6.8"
}
},
"nbformat": 4,
"nbformat_minor": 2
}