{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# 0. はじめに\nこのノートは「Open Problems - Multimodal Single-Cell Integration」の課題を理解するために、母国語である日本語で内容を整理する目的で作成しています。専門用語が多いので調べながら進めていきます。理解が間違っていることに気づいたりしたら、コメントで是非教えてください。ちょっとずつ更新しています。最終更新日は9月27日。\n<br>\n<br>\n<br>\n","metadata":{}},{"cell_type":"markdown","source":"## 課題要約\n* 細胞ではDNA→RNA→たんぱく質の順に情報が翻訳されて物質が生成されます\n* `train_multi_inputs.h5`というファイルを使ってDNAに関連する量からRNAに関連する量を予測するモデルを作ります\n* `train_cite_inputs.h5`というファイルを使ってRNAに関連する量からたんぱく質に関連する量を予測するモデルを作ります\n* 課題の意義としては、測定量を使って幹細胞が何の役割を持った細胞に分化するかしれたらいいなぁみたいな気持ち\n<br>\n<br>\n<br>","metadata":{}},{"cell_type":"markdown","source":"# 1. OverView\n\n## 1.1 Description\n### Goal of the Competition（コンペの目標）\nこのコンペの目的は、**骨髄幹細胞がより成熟した血液細胞へと成長する過程で、DNA、RNA、タンパク質の測定値が単一細胞内でどのように変化するかを予測する**ことです。競技参加者が行うことは、4人のヒトドナーから5つの時点で採取したCD34+造血幹・前駆細胞（HSPC）の30万細胞の時系列データセットで学習したモデルを開発することになります。このデータは創薬企業Cellarity社がこのコンペのために作成したものです。\n\n> **_NOTE:_**  骨髄幹細胞は造血幹細胞のことを指しているようです。造血幹細胞は分化すると最終的に赤血球や白血球、血小板などの血液関連の細胞になります。\n\n<br>\n<br>\n\n\nテスト用のデータセットはデータセット内の未知の後期時点から採取されています。競技参加者は1つの様式(modality)を与えられ、同じ細胞で測定された対になる様式を予測することが課されることになります。このコンペティションの挑戦は、テストデータがトレーニングデータのどの時点よりも遅い時点のものであることです。\n\n> **_NOTE:_**  様式(modality)はゲノムモダリティのことでしょうか。ゲノムモダリティとは、ゲノムの構造や機能の事のようです。図を見ると、測定量から遺伝子発現量やたんぱく質量を予測するという課題のようです。細胞の学習していない細胞後期の状態を予測するところがチャレンジングな点です。\n\n<br>\n<br>\n\nこの解析により、細胞の状態を何層にも分けて遺伝情報をマッピングする手法の革新を加速させると見込まれます。もしある様式を別の様式から予測することができれば、これらの複雑な制御プロセスを支配する規則についての理解を深めることができるかもしれません。\n<br>\n<br>\n\n\n![image.png](attachment:d9301e6f-7a9f-4b1a-83b6-fa50acd653e8.png)\n\n<br>\n<br>\n<br>\n\n### Context(背景)\nこの過去10年において、シングルセルゲノミクスの登場により単一細胞内のDNA、RNA、タンパク質の測定が可能になりました。そしてこれらの技術によってこれまでにないスケールと解像度で生物学を研究出来るなりました。その結果、ヒトの初期胚発生の詳細なマップ、新しい疾患関連細胞の発見、細胞を標的とした治療介入などが可能となりました。さらに、近年の実験技術の進歩により、同一細胞内で複数のゲノムモダリティを測定することが可能になりました。\n\n\nこのようにマルチモーダルな単一細胞データが利用できるようになった一方で、データ解析の手法はまだ乏しいのが現状です。単一細胞の体積は小さいために測定値はまばらでノイズが多く、細胞間の分子サンプリング深さの違い（シーケンシング深さ）やバッチで細胞を扱うことによる技術的影響（バッチ効果）は生物学的な差異を分からなくしてしまう事が多々あります。マルチモーダルデータを解析する場合、異なる特徴空間やモダリティ間およびバッチ間における共有および固有の変動を考慮する必要があります。さらに、現在の単一細胞データ解析経路はたとえその下にダイナミックな生物学的プロセスがあったとしても細胞を静的なスナップショットとして扱っています。時間の経過に伴う状態の変化と同時に時間的なダイナミクスを考慮することは、単一細胞のデータサイエンスにおける未解決の課題です。\n\n\n一般に、遺伝情報はDNA→RNA→たんぱく質と流れます。DNAはRNAを生成するために読み取り可能でなければならず（ATACデータ）、RNAはたんぱく質を生成するためにテンプレートとして使用されます（ADTデータ）。これらのプロセスはフィードバックによって制御されています。例えば、たんぱく質がDNAと結合することでより多くのRNAが生成されるのを防ぐことができます。このような遺伝子の制御が、生物が発展し環境変化に適応していくための動的な細胞プロセスの基盤となっています。単一細胞データ科学では、動的プロセスは生物学的プロセスの進行を捕らえるためのいわゆる擬似時間アルゴリズムによってモデル化されてきました。しかし、これらのアルゴリズムを一般化し擬似時間と実時間の両方を考慮することはまだ未解決の問題です。\n\n\nコンペティション・ホストのOpen Problems in Single-Cell Analysisは、単一細胞解析手法のベンチマークを標準化するためのオープンソースでコミュニティ主導の取り組みです。Open Problemsの中核となる取り組みには、既存の課題を測定可能なタスクに形式化すること、高品質のデータセットを収集すること、コミュニティから提供された手法のベンチマークを集中的に行うこと、多様な手法開発者を集めて単一細胞のアルゴリズムを改善するコミュニティに焦点を当てたイベントを開催することなどがあります。Cellarity、Chan Zuckerbeg Biohub、Chan Zuckerberg Initiative、Helmholtz Munich、Yaleと提携し、学際的なコラボレーションを通じて遺伝的動態の時間的変化の予測にどんな進展が見られるか期待してます。\n\n人体には約37兆個の細胞が存在しそれぞれが異なる行動や機能をもっています。一つのゲノムがどのようにして多様な細胞状態を生み出すのかを理解することは、健康や病気において組織がどのように機能し、どのように誤作動を起こすのか、メカニズムに対する洞察を得るための鍵となります。競技参加者はこのコンペに参加することによりこの単一細胞生物学の基本的な課題の解決に貢献することができます。予測問題を長期にわたって解くことができれば、血液や免疫細胞の成熟に伴う分化に遺伝子制御がどのように影響するかについて新たな知見を得ることができるかもしれません。\n\nコンペティションのヘッダー画像：Pawel Czerwinski on Unsplash\n\n\n> **_NOTE:_**  単一細胞測定が発達したものの、測定ノイズが各細胞間の違いよりも大きいことがあり解析を難しくしています。また、細胞は常に変化しているため、ある時刻の状態よりもダイナミクスのほうがより情報があり、個体スケールで現れる適応の基盤にもなっています。細胞ダイナミクスを実時間でもモデル化することで理解できることが多くなるようです（具体的にはよくわからないですが）。\n\n<br>\n<br>\n<br>","metadata":{},"attachments":{"d9301e6f-7a9f-4b1a-83b6-fa50acd653e8.png":{"image/png":"iVBORw0KGgoAAAANSUhEUgAAA+kAAAH1CAYAAACKtr+2AAAgAElEQVR4nOzde3Qc53nn+V91gzdbFglaiZyRZJKQPB5fIsKB4sQnckhlwPGcaBIxcSB7nDNHnKwWXF9kc042AZO1IyrxRMQkWdOWHYW0osuOJdtkZEvjyCuH0CEV22s7JixKkeibiAaakCVRJBq8E7eu/eOtV11dqO6u7q7u6ga+n3P6AF1dVf2CBKrred/nfV4JAAAAAAAAAAAAAAAAWGj6JbkxPQaa3PYojqq4jfuTbQ4AAADQGKmkGwAAAAAAFexXfQMME2WO721Ae8N0Stob4/l2q/jnOBrjuZEggnQAAAAA7axSkN0nEyAnqVfSIa8tQFkE6QAAAADaWa/KB+E9zWpICbtlMgG6Em4HAKBF9Gp+Wld/oi2qHnPSAQBY3Mqlu7sqP0J9qMKxjUx3D6sbBJTFSDoAAACAdldqtLyzzGtASyJIB1BKr0whllKFVnYq2ryqPm/fsHMc8t4jjvSvvhLvwQczAAAL04jv+1Kj4cHtww1qCwAAkdWS7h6sFlruUSr1vFOVU9OCQX8pldLdO0P2qZT6BgAA2kfYPcXewPOweenBgYKwe5xK6e69Ieex5yp1rzEQsn+lVPvgvUy/zGCD/2efUKGafbXV3Tu9Y8PS/3er/aZDAkDbqjZILzXqXe6xO+Q81QTo9lFqCZVKQXrwQ7rcuQAAQPsJu68IBsJhAbM/ID2k6oP0sHuM4OOQ5mcFxhGkl8ponPD2ryZID5sbX+pnIQsRABqsmiC9K2TfYDXSXs3/EJkoPk2k9+zT/A+eCYX3gpcL0sM+BMuNygMAgPYTFqQH71uCa5D3aP79QTVBejUDDkdVfL8UR5Be6mHvc6IG6VHb4r8foxI9ADRQtSPpXTIX890qHTSHXez9va5hqWXl2jbgPUotn1IqSA9++JZ7LwAA0L7CAmapeKQ8OGgQNtIedg8TFqSH7efP0utTtNVnqqnuHjVIt/dcUYL0sPvA/b5zdCo8W4D7KQBooEYswRb2geP/gAt+aByq8/3CPgTD5qGzNBsAAAtTqSA9ODDgHzQIm7MeJUjv1Pxsv7B7p7AMxOC56g3SJ1R63nuUID3473ZU4YMipTIVkACquwOoRr8Ko+zljASe96hQ6CSutUh3q/jDIyfp5pjODQAA2kPwnqO3xPfDMvcKUfSrOJDNSdpT4r33lXn/OAyGvEdUXZrfnj0K/3cY8rZv9x6rNf/fFk3SkXQDALQ0m4JebRG2fZo/L7wzsG3Ie+xT9R8CYR+Adh3UoSrPBQAA2ldwSbUe31d/oF3N/UGwcFq5ZduGVTzSHXeQHtY5EFVYW0r9Owx6D7QARtIBhLHpUztVW5X0EUlbK+xjlzM5KpMOH8dyaVR0BwBgcbCjwcER8t7AV/n2iyqY5h02ddA+goMScVZGH1H00f8wYWntjI63AYJ0AEGHVHrO+qBMClSlAFwyPb+bFO1DsUdm3li9RUp6xdroAAAsBv7g1T863ClzP1AqSI8S9JYqZBtVvcdb9QToUvic8nrPiSYgSAfg16/5PcBDkq6T5MgE6NWkQtljr/aOrRSw96u6onabNL9HmOXXAABYXIL3F30qDtJzau4IcqsE6WhTBOkA/IK9ziOKPhpezohMcG+D/a0yQXvYh0/UuVxbZToBgp0GXaq/ej0AAGgfwXnWwfuAeuvV7JG5f4n6aJWU8rB2ULG9DRCkA/ALXrhLfcjUe4HfIxNcr9b8QD1K7/OQCoVU9ih8ND2uXmwAANDaKlVur3awoVQxunYT9m9CkN4GCNIB+AUv5mEfSp2KNlI94D2Oqrq55rWkdgXnyEdtIwAAaE/BDvpyo9f+kfQo9xlhQXo7dv6HZRCU6nCwa77bosFxV6kHAPiEVSQtFcDuDNl3r+/1PpmgO6y6qS3Y1imzJnpY9VP/B1yXd+7gfmEfCsH33B+yz/7APhOitxgAgIUi+DkfvBcIu4ex9wN+YfdFwXuPrpB9wmre2HuevTKBbVjx2v6Qc5W6P4lyv+O3W8X7Hw3ZJ+z+KKzDYSCknayaAwANUk2Q3qXwADvKw38h76vxHP4OAb8oH1phP2e91eIBAEBrqBSkl7r3CN5bRAnSpfkBsL2vsAFuV0ibwgLbsPcb8J3DH7A3IkgPe/9DgZ85LEAvFcwDAGJQTZBeav+wD6lgj3XYh2U1AX+5quxRP7TCRuYZTQcAoP1VCtI7FX5/ESVoDgvSO2WC2WoGGw6VOE+5Y/z3ZI0I0qXwILzSgyVtE8ScdABBQzJLpoUttbZdpkL7Vs2fr9Wr4oB4n0xhuO3eI8yI75yl9qlG2DlYkg0AgIUvp/ACcbVWds/JrHAT9fghb/+w85RbvrYZo9WDml+/p5ScpJtl7uMAAAAAAAgVHEkPm9IWNi89GAT3hOxTqUhar0qnhEctstav+SPz+wPHNmok3er02huWIWDn1QMAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAChBiRNSHK9x4S3rbPG89nz9HrPe33bqmmTK2l/jW0AsHD1StqrwnXF9Z731Xi+/d45BipsK6dT0lEVX/sAAACAqnSpcFMZ9jjq7VOteoN0//4E6QD8bAdeqUfUoNovjiDd7k+QDgAAWkIq6QagJrtlgvARSZskOd5jj/d6l2q74a3HThGYAwjXK3ONkMx1yl6zNslcx+S9XkvnYj1tOiQCcwAAANTJP1oddnO5XyZAD97s9qs4NX635qfF1zqSvluFEfxDYiQdQDE7Wh12Xej1tvcHtnepODX+qOanxdc6kt7lO6//PQjYAQAAUDWbMnq0imN2qnRavD9QrydIH/B9T5AOwM9eR4KBeCldKu5ULJUWX0+Qvl/mGucP2AnSAQBA4kh3bz92hHyk7F7F+9ubVZsav1rSsPda1JvmcrZKGozhPAAWHn9WT9Trli2AOaRCavzWwGv1sFOFhuo8DwAAQOwI0hc+/8iQHWGakNQT8joAtIKwbJ7d3rZOFa5fAAAACw5BevvJeV+jFliqNOLUzEJNABafnO/7qNebSvtx3QIAAAsWQXr7semZXQofBd8tMwfdjjT5b5CdkMfVjWkmAEgy1yB73QpbD71LpuCkfw65TYvfrvDr1h4BAAAsUATp7WdIZj65ZAJyf6A+IDPHfECF5Y78cy7tNn9RpmYv1QZg8bFBda8KaeuS6Uzc633dqcL1zF63+lXocLQFMCca2lIAAACgBl0yldnDKh/bqu3+dNBGV3f3o7o7gDB2ZYpSD3+HYaOru/tR3R0AALQURtLb04hMmvp2Faez57xt16m4ivL2kH33yFQ39m8DgEYZlLnm7Ats3yfpZhWvEDESsu+IWEkCAAAAAAAAAAA0i5N0AwAAAIAYdUta5X0/6j0AAAAAAEACDqpQa2JHwm0BgKoxJx0AAAAAgBZBkA4AAAAAQIsgSAcAAAAAoEUQpAMAAAAA0CII0gEAAAAAaBEE6QAAAAAAtAiCdAAAAAAAWgRBOgAAAAAALYIgHQAAAADQitwWeuxo8M/6qo5mvREAAAAAAFVykm6AmhigS4ykAwAAAADQMgjSAQAAAABoEQTpAAAAAAC0CIJ0AAAAAABaBEE6AAAAAAAtgiAdAAAAAIAWQZAOAAAAAECLIEgHAAAAAAAAAABA7A5Kcr3HjoTbAqA+bkLvO6rCdaTU43Cj3pyRdAAAAAAACqJ08G1r1JsTpAMAAAAAUHC/pCfLvP6oTNZOQxCkAwAAAABQrNxIecNG0SWCdAAAAAAAgg5LeiBk+x0yc9YbhiC9yXIPPrZm4u/3Pp50OwAAAAAAZW2TdMr3/JSkXY1+U4L0JsvPXfyd1MXpNybdDgAAAABAWZMqDsq3edsaiiA9AalT538h6TYAAAAAACraIWlMppDc/c14w45mvAkKUudnL0m6DQAAAACAyJoygm4RpDdbfq436SYAAAAAACJ7pJlvRrp7AlIXLl6adBsAAAAAAK2HID0BqTPn+XcHAAAAAMxDsAgAAAAAQIsgSG+y1JnzXUm3AQAAAADQmigc12TOXP51SbcBAFrR3d98/8ZXn3Ro9IPv+tJorefa8+2+7rl8epV9/sF3f+lgnc0DAAALxCefyt7ruvqlRr7HXN79zo7r1nywlmOTCNIPStpYca8FbuKh/3X96g/89reSbgeAROySWWfzcNINaSWOkz9gv1+7ev3YsbO3FQXpjx3+9IZSx97Y/bEn/c+//ZMvrp88f3yVb5MTW0MBAEBbW+o4v/GLncvWXLYs3ZDzn5ia06GJizUfn0SQXvImazFxZ2euk0SQDixO3ZJWVdxrEfrVN/2uXn/JlZK0RnLX+F+78R0fLXOkW/TZ8mtvfr8k6eTZcX33p1+Ju5kAAKDNXbYsrStWLEm6GaGYkw4AAAAAQIsgSG+yjuMTjJ4BQIhrr7rhgUtX/Fys57x0xc/p2qtueCDWkwIAADQQheMAAC3hqsvePhr31PEl6WXeeQEAANoDI+kAAAAAALQIgvSEpC7MUEAPAAAAAFCEIL2JJv72oWtffTKX70ywKQDQch576q7b2+m8AAAAjcCc9CZyL1n2G0m3oVVls9n7Jd0i6VQqldp45ZVXsn40AAAAgEWHIB2Jy2azW2QCdElamc/nD2YymbXr1q2bTLJdAAAgkk5Jfd7XqIYlDTWmOQDQ3gjS0Qq2BZ6vTKVS2yTtSKIxAAAgsgFJO2s4brsI0gEgFHPSm8hdsvSSpNvQajKZzFpJ64PbHccJBu4AAKC11BqgAwDKYCS9idK5M73NfL9Hh3dcOyf3Y66bf9usppantWx0WXrFnt96x59+vZntKCeVSm0u8dLKbDa7+Y1vfOMjTW0QgMTc+I7b7pAUe5E3c969cZ8WgNTv+35I0j5JIxGPjbofACw6bRWkj4+Pd+fz+YOS1O7FxdK5M29r5Pm/OnzH4y9eeO49gc3rJd304Pc++nJKy/r+86/81Tcb2YYoHMfZWOblzZII0gEAaD09krq874ckbUqwLQCwoLRNunsmk1nlBegrZYqL7Uq6TXVxnIZ1kPzDoY8///L8AP1VZ2ZeuvzUzNg/7/2XgT9rVBuqUClIBwAArcdfJG5PYq0AgAWobYJ0r5DYSt+mDceOHSsX4LUcN51a3uj3+IdDf3rkxMWfXB1l34npzB2PDu+4tvKejTE+Pt6t4v/ToJXePgAWAdZJB9pKrsT3AIA6tUW6eyaTWRVWSMx13S2SDibQpJp0nDz95kae/8Cz997648nH31LNMTPu3J8roRHrfD5fMQCfm5vbKKltpzUAwMe/n9noyjlgn19zydIn/8MbXlv02fW3z+dKdiR86JrOO/zPHxw7veXUzNwaSVri6Hs7rlv7q3G3GYhgWCY471Qh7R0AEIO2CNLT6fRmhY+43pLJZLaxnrZx6ZIrtld7TG5q9D81oi1RuK671nGcsvt4c9bbe2oDAEjafOXrJEmXdqQ3OI6zwf/ah9+0utyhRQH8f/yFSzSVz+tHp6c1em6m4RlaQBmDMtXdB0TKOwDEpl3S3beUesEL4Be98fHx7pSzNFKau9+cO5X+0r/8YX/lPeNXoWicleiUhvueGf3Nz//r2ON//69jn3n4yMiaJNsCoL1dsWKJrlixRK9bUt9H72XL0uY8He3yEY4FbFAmOO+SWUKBEXUAiEHLj6R762hvKLPLRkn3N6k5dXHm5l7TqHPPzc3V3FnhSrcqgR5wx3HWuq5babeVmUxm7bp160ab0Sa/zz899tTzU2635ErSe5Y6qY/c++zon/3B29d+stltAdC+Zmc12n3Z8jFJsXb0XfGaDh2fmn08znMCVdrv+77Pe4wo2vJq+8ToOwCEavkgvcw62lbbjKSnT55a2qhzO45T87/D1Nyp9XG2JSrXdSPdsKbT6W5JTQ3S//5fRw+OXHSL5sxPu3KOXtBf3PvsqAjUAUS1813rRn82eXpUcQfpK5boyuVLCdKRpN6QbV2KNqI+FHNb2lUr1Va6X20y8AUsdC0fpEdIiV45Pj7e3W5rpnccn1gV17kymcwqmTXQa3Jx7vTSL37vj97dzHXTjx07tjHCKLokyXXdbjVxvfS9z41d//R5t2T2xrGL+vOHj4z8z/e+tWusWW0CFoMb33HbHQrMv47vvHvjPi0AxGGDpBuSboTM1NK1STcCgNHyQbqkmyrt4KV6t1WQHqeOjo7uqAFvKcs6XrNFUtOCdNd1I3dSRJy7HpuTs+V7kaddOcdn0k9IuqZJTQIAoBWVr/6KqFphNL2tljUGFrqWrjqTzWYjpXA3O4irxcn/5+F5S8jFJZ/P1/3zT+fPb4qjLVF5o+NRVbVW+gPPjt76fz81+srf/+to1R96Dx8ZWTM+7VYswDc+7V6997mx66s9P4DF6XM/nShXW6Vmn33+5IHKewEAgHbS0kF6FcH3Bi/le1GKo5Pi7MxLV8TRlqgcx6kmpWrl+Ph4pEB975Hs+0Yvas/JWV02clEb/u7p0e9W064zc6nIy72dmtVd1ZwbQHmPPXVX7KnujTwvAABAI7R0kO66buTgM51Ot/xoetDE3335N2M6VVUjzWFm81OpL37vj94dR2Miqmre09zcXKT/3xem83837RbS745N61eqGfF+ZVbvibrvz2bcRAruAQDQwjplCsrZB8uyAUCVWjZI95ZeqyYIaukq787M3H8MbnNf0/Fv6z2v9++0st7zSK/OS2+WqjoWomQL7D2Sfd+JGc3LqKg0x9x6+MjImslZrYjapmlXzr3Pjn486v4AACxgA5IOSZqQWZrNPo562/YqvBo8ACCgJYL0bDa7eWxsbIc/Zb3akXF/EJfJZFZls9nD2Wy2ZYrJOVOzy4PbUudnL6n3vB0dHbFV4pyaO9vMbIRqOxYqtu3kjPvfw7Yfn3G7djxVeTrEebfjE1W2SXNzznurPQYAgAWkUyY43ympp8w+fTJB+0CT2gUAbSvxIN0bCf6q4zi3p9Ppg75AvaqA0XXdNceOHdsoSel0+hGZUfj12Wy2maPDJbnp1LwgXfm5unuU4ygaZ52dOR7rGr6lRJ1fHrDS/v+WcnzGDU2pm3blXL00/deV3uDsbPTpFdZk3n1ztccAALBAdMoE3v7gfEjSdt9jn6Sc7/Wdkvqb1UAseP0yv4Ou9ziq1ukIOqrq1//sV3H7d3vnwSKTeJCeTqd3+J6u7+jo2CHVVgzNdd1t2Wz2fpk1J+WdZ0eZQ5qm4+TphgRzjuPUPR/dmnOn0l976i/jmidfUjXLr/mV65D4wrNjd/rnogedms3fWOn8L8+466pt0+SsVkQZpQewuH34TaufbMR5P3LN61thfWUsXv0qBOjDkq6WtEnSoO9xs6TVMgG7tVMmwF9sNkraEXgoZFvb1VlKyG6Z36UhmeUAHZnfuZ1qjUC9U+bvIqoumZ/J36m1VebvCotMouuke6Pot/i3ua77sWw2u9Z13VpGdeetqe667prx8fHuK6+8smVS363UmfN1F1NxHKfuNdL9puYu9Ev6emwnDJHP5zc6TvVLq5bruDk1N///3u/Fab3h4SMja9771q6xsNf3Hsm+7+lz+Zo6rf7tstQWSZGrwgMId+M7brtDUuyV2M15qx3MABCBDYRGZILzXJl9B72vNkDvk7SncU1rSaOStmn+lL/gde+R5jSnrQ3I/A5tUnEgvEcm2B3wvi/3O9lIvao+SLcZttUcgwUq0ZH0VCpVKhW9bMBVLdd1E095T586M2+01ZnLv66ec2YymVU1dmaUNDV36lfjPF/MSi61d2rWfVOlg8vNOT8/l39/rY06N6d5RQEBAFjgelQYDY8aDA369luMReRGVblT/9OSWm5gqcV0qhCEhwW0e1ToCLL6ZdLGbVr87sAxXd72fpkaC67M76hNn7cp9bZjqs+3nyvT+eRnM0xGfNsGZIoo2mOCae22TYdUyFKxbar25+iT6Z12vffsE9pKqwbpsapmKbdGcaZm5m1LT565tJ5zdnR0xJbqbp2Zeeny+w5sbmj6dj3ruocVFNx7JPu+c/nKWSHlUt4nZhV5mbag2byuqfVYAIvD5346saHyXtX77PMnDzTivEAE/nT1akb+7L6LMd1dMunsoVl9kk6pkAKP0vplfn9KZWKMyKS+29cHZILZQW/71TIBuD/FygbV/TIp5o7M76odld+nQjp9v3es3Xadt80fqHd57bBB+oD32OQ7xm6T955DKqTu7/G1aaiKn8N2fvX59hvS/E4EtLjEgvRsNrs57lHgMtZ7qfWJOP7IP4XO83amZ+r698/n87EH6ZL02te++aONOK9PzZ0AYQH+qRl3e9i+QTblPbj94SMja07O6rJa23R2zr281mMBFDz21F2xp7o38rwAXlVNwL1Yg3O/bSW275A02cyGtKkeFQfA5dhR9+0qBO0j3vd9Kvw+2imogyp0JNkgeY/v2E6ZgHefClM4hlUI3q1eFXdIDQTOPSyTVeIvutijQkBun+e89kb9Oewx233vNeL9fPzttZEkR9KbmoJe7ZJusb736TMl10Of3P2lW2s9bz0j0uXMafaWynvVZX0dx24Obnhl1n171IPDUt5rWXrNb3JOdS+lBwBAm/GPnkdNXe9ScaG5xeoRScFikmOivk1UXYr++2MD5+Couw3wbXDe623b59vHH6T7z1dqFL9ThUDYPx89J1M80Qb1fSrUZvC3o1PF00b8P2c1P8eQijswbKdGUvPzUYNEgnRvVDvWeecRJBakOzNzJecsu3O6tubzxljZ3e/U9LGualLe7/7O+yNnKdSb0WALAdrnX3h27M4oqe7WsYv5eR0QJ2byda91vvdI9n31ngMAgDaSU2HUr1+Vl1XrVHFa7lCpHReJ4Gh64vWTFqhSAaoNpu32Ls0PgLsUHvBKxcu+BeekB4vG9agwf90uWWgD+hHfMdL8kXT/OSr9HD2+NvsFR+jRBhIJ0ps1F92vUaPOkczmywS8+ZqC9EYUjfOLkvJ+34HNq774L3941JnLZz7/7f8y+8Xv/dG7Kx3T0dFR97SDfD7/6nytF6fdj1Rz7Lm8Ov7nc9l77POHjoxtOzFTe/q9NeVqU73nAACgzQz6vreFr4Kj6l0qFOPyz7Fd7EHDYZkicZL0qKSDCbal3QyrOE08aK9MMCwV5oYH2ZFzmwrepfkBsN0n7P2dEg9/Crs9dq+33ZGZkx5MRbftzGn+yLr/eaWfIyzQr6XKPFrAognSXdddU6oyeKOVWyM9PXm2ptTvRqfvz7gXPlRpn9SKy4+emj7WJZk11vOa2lfpmHJrnVfhpm/8eOz+zz49+mItqeYjF/J/8MCzo7fe98zob/7kvPs3MbRHZ2eTL04IAECTDal4/fN+zR9hPCoTvNt03JzM2ukwc9BPqfQcdYQbkvl9Cptm0atC0TQpPKC3x+7zPbfn9QsLbnOaP7c7WIXdXzSu1Oh2v4qzUfyj5vbn8Lcpys/hD9j9bfOfB20ikSC9iQXjitQT2OYefGzNxL1f+dGpv77XnbzrCy+VKgYXJnXhYskq7h3HJ2rtOJg3NztOZ2ZeurzcyPg/HPr48xdmc6uDx3xl+M/+W7nzOo4TSwG/c3O65cVpvaGWY6ddOT+5oM8/P6XHZtx4/gZys8n8TgNoDx9+0+rg/NNYfOSa19/QiPMCVRiUCbqjzHfdJ1ORmrmxxqSkbpml2RDdPpmgc7eKlxazS6VtVyEo3adChXZ53++VCWRtB5O/QJtlR6CDo9f+ddjtsXu999vjO9ZfIM4/yi1v/7CAOzgf3d+mqD9HMBgPC9zRBhJdgq3Zak15z2QyqzpGsj9Z/sxP3rxk/LiWHh2//DXDP/palJH5TCazKnXmfNl/51qKxzUjfX9OF74Ytv2rw3c8fuLiT64Oe+3kxcxfPXHknpIBa6LTDhroXF4dDzw7WnMRQCCqc9vu3GEfF7bdWVdW0oVtd270ny+uNtbqxnfcdkc7nRfAq/bJFMa6WYUAyT72eNuuVvRgfjEhQK/NJpnfrd0qZG30yvyO+adhDMksb9avQmbHsMwSaFZYITobJAeD3mHvPez5Dnnt8E97DI7A36zC+uWud97g/Pd9Mh0OE97z4Mh6pZ+jVFo789HbVOSCWwtBreulr3zi+wfTJ08t9W9LnTmfWjU0/NeSygZml37zB5Vvos9e+C+S7qm4n2d8fLw7n883fOT27MzxK/b+y8Cf3fzOwT+32/7x8N/cM372e+8pdcycO5V+8czhb0u6Mvjak0ce/M2ksiia4fycs11V/D8CtXDk3K7LOqVlS+Re1qnzR7P3Fe1w10OlD77tA0VP3R+OyPnhiDQ1I53ISazPC6A+Nu12sOxeQDwGFe13zb+EWpitVR6zT8VV4INWB54Py3RSlXvf7SqeOhJW66hcm+z66pXagjaxqIJ0mfXSV61bty7yGpQvD/3ztcv2fz903njHS6/8F1UI0p28885K79GRO/2rUdsjSa7rNm1O/8R05o5/OPSJX7pkyc8dzE2N3TZ+9ntdlY45O3P8in849IlHfu+6vyhKyX/l/A8/s+6SirXl2tb4tHv1A8+O3nrL29cSqKOxfv2XpCsuD38tEIiX9ZYu83jhZekrT8TTNoT63E8nNnz4TfHfK332+ZMHFH5jBgAA2tSiSneXpHQ6XdVc7iUv5f7cmZoJP9fJU0tPfuGR8ktvzc6+tdJ7pM5eWHrs2LFIo/xeVfemFt47cfHHN42e+TOufpwAACAASURBVNanbJG4qMd85dCOV5ej+MfDf3NPqRT5heSFKd398JGRBZstADTSY0/ddXs7nRdYRHbKzPXdH9i+v45HpSXbAGDRWmwj6ZK0K5PJPBJlND2TyazqePAb/6ncPunc2e2Svlzq9Y7jk5ECNtd1D4yNjd2RTqcfufLKKw+X2i+VSm2TtDLKOZN2/OKRgQe+u/V/l6Txs99bFOk25/LqGJlK/fDhIyNvee9bu8aSbg8WHveGdz7gXNZ5S6wnvaxT7g3vfEC7Yj0rgIWjR6UradeKebIAUMJiDNJXptPpg5lMZvO6devKFuu45Ps//MeOE5PpcvssOfZy2SXUqqne7jjO7fl8/vZsNitJthLwqOu6o47jjLquu9ZxnLYaEQpWgF8MJme14tm5VObeZ0f/7A/evvaTSbcHC8zbrx6NPbt52VLvvAAAAEjaYgzSJWl9Op0+PDY2tmPNmjXzxo6y2ezmFc+NfPK13/ju2yqdyJmecSbu+8qdq//r7/5J8LWT+x6/Vd8uOSheyQb71XHMDbn9itY37co5ekF/8anDo3/437rXBtfTBGrm3PXF26uad17NeSkcByBcWBEriXoI81zYducWV2rotMS80tsu2fXHNd9gAmh9izVIl6SVjuN86tixY9vy+fz9juMc9pYH2+y67hrn9LnIJ0qfmNwiaV6Q7kxMsiTXIndiRpEzKQAAQHvLS2sdORsq71m7tOa4twAWuEVXOC7Idd01Xgr5V13X/VgtS4QteeH4GyYe+l/X+7flHnxszdLx478cW0MBYIFbyOukf/hNq5+svFf1PnLN629oxHkBAEByFn2QHpf0y7n77fe5Bx9b0/H82A9TZ87z7wsAABaDLknB6V29kg5JciVNSNrt7QcAKIMgMiZLxl68evKuLzw1sftL313yo6Mj6dzpFUm3CcDC4972nxsyKtyo8wJY8HpkllQ76n1v9Xrb7bZOmWXXDgX2AwAEEKTHaOnR8e7lPxz9FUbQATSK8/Vvb1QVNTMiOX3OnDdhPxj9xsYL06djPeeF6dP6weg3Ev3Ztn8ns/afXzm/Nu7zvnBhRl978fR/jPu8QBW6ZAJxuxSbP/ge8H0/ImnY+75T0t7GNw0A2hfBJAC0k6PZDTpzNt5znjlrzpuwF3M/2XA+5iD9/PRpvZj7SaI/W0eH1v7r5MWq651U8sL5Wb18cY4gHUkaUCHF3R+Id6kQuOckXec9tvte72tSGwGg7SRS3f073/lOEm9bleVHMlo2PtK09zsVw7/J6emf6SeTx2NoTeN0vhDP//2PL0hjF91YztUEwVG8SUksnQIAaHc20B6WWaYtF9guSXt82wdVCOx7JO1rQhsBoO00PUhfvXr10wMDAy2/dMRy11m+LK9lcZ8377j5GceZCW6/eOifLtZ77lTK7XCWzLy23vM00lemBk7FcqIlS5e7HUtj//9pkODa04clbUuiIWgpB2o56DW7tku74m7KqxLt+fo/rt+rBmbBJvazffKX10mS/ntjTr9e5X+2OzT/GgTEoUeFUXR/IC4VB+nDKjYsM8rOvHQAKKHpQfrExET3xMREs98WAFpJzXOkL2y789Vj55SevGTXH9eclXFh251rJb06V3rFrj85WOu54nD3N99f+Hfp0OgH3/Wl0VrPtefbfd1z+fSrHcIffPeXEvvZdjyVWTU7q277fHZWozvfta7mn+1P/yXTnXK0SpLyrib/8p3ryMxBEvyV3EcC2/0B+FBzmgMAC0ci6e4AgNrEGUiv2PUno5JqDhbjFmcg3f9r+1omcN3xjnWTkmL72QjK0eJ6fd8Pq3iEXZq/TBsAIIDCcQAAAKiWf/S8p8T3wTnnXb7XGWEHgBII0gEAAFCtERUC9X4VAvB+3z7+QLxLxUUnmledFwDaDEE6AAAAarHH+9ol6aikQyqksw+pUDSu33vdP4pOZXcAKIEgHQAAALUYVCFQ98upsCa6VDwPfVjSzY1sFAC0O4J0AFh4emWW5dqZdENiMiHz8wQfEzJrLi8EAzKjkP6fbaH8/2Fh2+o99smMkA9KulrFS68Ne4/tKl5PHQAQgiC9ftwMtx9uhrHQ2TmhfWX3an+dMn+7u5NuSB26ZK5HO1VccKtThWsV1bDR6vbIjI5vkgnEg0H4kKTrZAJ4AnQAqIAgvX7cDLcPboaxGHTKdB4OyfzOL5TOtWFJTuBxnUzxKVu0qh3tlbke2SDG/7MNea+183UXAABUiSC9PtwMtxduhrEY9MtcmwZl/mYXcgfisArFp9qxg61fhWvSJs1PD7bb+tS+110sHl2a/3fYq0Lm2oTMZyy/ywBQAUF6fbgZbh/cDGOx6Je5HtnqyT0yN8oLkf3ZRlT8N90ubEbPYJl9tsp0JrJcFVpVj6T9Kq7eLpm/zf2+bZ0y16dDgf0AAAEE6fXhZrh9cDOMxcB2Mtlqy7ZjrT9897bSo/m1Mg7J/LztWim6R2Z+7lCZfdrxeovFo0smELf3Pv7g259d6L936FTxeukAgACC9NpxM9xeuBnGYtAv83tur0vDMr/zvWrPDJgo7LQjAM03oMK1xR+Id6nwd5mTmVZ2nQrLsnVpYWcfAkBdCNJrx80wgFZis106VbxKg93W7h2IYbUyVstcd3eKG34gCfbvbliF+i7+7ZK5T7IV3f3V3Ul5B4ASCNJrw80wgFZT6e+y3a9LYXIy01Sk9rwuDSt652c7/nxY2HpUGJTwB+JS8e9rMEvNPidIB4ASCNJrw81w++FmGAtdv8I72ByZ0SvSS1uPHXUstzLIgEwnMNBq/FmDI4Ht/gC83DQzAEAIgvTacDPcfrgZxkJm54XuKfH6QqqZ4depwrKJ7VhTYp9Mu20VbP/nhq2YvVPFq2sArc7fGT6s4hF2aeFOCQSA2HQk3YA2FOVmeEDmZngh3VQttJvhPSr8//TI3Aj3ipthtCe70kSp65L9ve6T+X1vx79hW9AyTLmfvdXdrEJ17LBMnxG1b8FOLGz+0XO7xKn93gp+nnb5XmeEPVmPSFrle77W+7rRt21S0uamtQjAqxhJr17Um+Fete98q7Dq7hMqLMHWzjfDIzI/x14VV663Pxs3w2g3wZUmSrE3ywsty2ePpE2aP1rXLkYkXa3iglrWoEwxLpaFRCsaUeF3s1+FANyfseMPxLtUvPQav9fJOixpg++xxnv4tx1OrHXAIkeQXh1uhrkZBlrNPhWm2kTZb3uF/VrRaoVPL3JkamUshL/b7Zr/c25X+15vsTjY+6EuSUdlOr1tOvuQClk7/d7r/lF0staStUvSWJnXxyTtaFJbAASQ7l4de5Mb136taHXSDWiC7WrPQAUAgFZi6/AE613kVPw565+HPiyy1lrBpKRtkr5a4vUtTWwLgACCdAAAANRqqwo1Xzq974MZa8MqTAcMLteG5Dwi6UmZ1Ha/JyUdbH5zAFgE6QAAAKjHHpWfCjgkM6UMrWebpKcC2xhFBxLGnHQAAADEoVNmRH3Ae/hXLAhbvQDJOyzpAd/zT0saTagtADwE6QAAAKiHXaZ1QmZJwZ3ew7/KzX7v0a4r3yxk2ySd8h4UiwNaAOnuAAAAqFWPTPDdWWafXt/X/TIrxQyX3h1NNilT7X3SewBIGCPpAAAAqEWnigP0PQqv3G6Lxtljdje+aajSDplAHUALIEgHAABALXaqEKBv9R5h65/nZIJ3+1qPpL6Gtw4A2hRBOgAAAGphA+0hla/ubvnXTmduOgCUQJAOAACAavWoMIoeNnoeZkSFuegE6QBQAkE6AAAAquUvFDdSxXG5uBsCAAsNQToAAACapVwVeACACNIBAABQPf8Sar0l9yrWpUKa+1C8zQGAhYN10gEAAFCtnEyg3SupX9KgKqey+5deY530Brqw7c4t7orl2xr6JrP52df81R9e19D3ABYpgnQAAADUYlAmSLfrpd+s8PnpPTLLtdkR92Exkt5QeWmts/aK9XrLusa9yVeeaNy5gUWOIB0AAAC1GJIJ1AdkAvGjKh4h7/ceXb5tOZn11NFor3utdMXlSbcCQA2Ykw4AAIBabVfp9c+7VBygj0jaJFLdAaAsgnQAAADUY1DS1TLBejCNPSezjvpWbx8CdACogHR3AAAA1KJThWJxIzLB+mByzQGAhYGRdAAAANTikMw89IGkGwIACwkj6QAAAKhWjwrzzSstvQYAqAIj6QAAAKhWp+/7sGXXAAA1IkgHAABAtYZVGEHvKbcjAKA6BOkAAACoVk5mObWcpJ0y89K7yh4BAIiEOekAAACoRa/M8mr9MoH6TpnU9yjp7/sk7Wlc0wCgfRGkAwAAoBY7Q7Z1KdqIenA9dQCAh3R3AAAAAABaBCPpAAAAqIWTdAMAYCFiJB0AAAAAgBZBkA4AAAAAQIsg3R0AAABx6JRZM92/bvqIzJrqUSq+AwBEkA4AAID69Misk95XZp9hSYMyS68BSdslqTvpRnju9x7AqwjSAQAAUKt+Sbsj7Ncjaa/M2uhbG9oioLJuSY9IOpxwO7ZIWptwG9CCCNIBAABQi14VB+gjMkH4sG9bj8y66f3e835JOUnbm9FAoIzDkg4m3IaNCb8/WhRBOgAAAGox4Pt+UOGB95D3dY+k/TLz1ge858xTB4AQVHcHAABAtbpkRtIlM8+80sj4sKSbfc/LzV8HgEWNIB0AAADV6vJ9H7UY3JAKqfC95XYEgMWMIB0AAAD1yDVoXwBYlAjSAQAAUK0hFQLuakbF7Rrqw2X3AoBFjCC9/Q1ImpDkeo8Jb1tnDefq9Z3Hf35XpthLtW3Zq8KHMQAAWFgGva/9ivZ5778/2dOQFgHAAkCQ3r66JB2VtFPFAXmnt+2QiueLNdrOkLb0yQT3zWwHgNbVK9N556q4M6/WAlJhnYj7vW0DoUcU2DWbgx2cAKIb9B6dMn97pQYJulS4T5DMOulUdgeAEgjS29dumQ+9EUmbJDnew/ZMd6l5N5x2ORXJfFg7kq6TSYPrFBVcAZhrxH7Nvx70yQTLzQyQu0LaYjs4d4YeASDMfpkOL/t5v1Omw+uQ99p+mQGFoyr8jedU6MQPe/QLABa5SkG6DfRqSZ1G4/SqMP9rqwprkPqfb1chDc3qV3E6+m7F83/rn4tm33NYhflmjKQDi1uvCsHvHhU6FTepMJq2U827VvTLXPtykq5WcQcnAQIQnb0fCd5L9PheC/5dd/peC3twzwlg0asUpPfI3DgdlRnp4OalNdh5XyMqDtCtTTLBsj+VbKfmB+X9Mr3d9X4g7lPhptsWkekMtBPA4mVH0IZkOhIVeG6/+q8VXSpORz+q+LJybMfiHt972iWkOkXHIgAASFBHhdftaGiPzM1Rn0ywN+Q9KPqRDHsDGTX49ae+b5L5v7Pzx3pkgvXgqHu99qowUsXvCdpRpc4rlhGKzgbFYWsp288Tvy7N70C0QXtYllC1rgvZZjsVc6JjEYjKSboBALAQVRpJH5G5mblO5sZoWIU5xrtlUqcZYW99/nR0W1RpQoWb0mqWToliv4rT8Qlm0E5s9tBEhcehpBrYZvyj0lGDXzvNakiFLJ2tgdfi5E/Hj7vDEgAAoCpRC8cNy9y4RAnYKRLWeDbojZqSWemGNq7UTjs67w/Qw0bOgFa1WyYIJN05Wb2+r/4aGlLxVJq43stWh7efdQCiqedvkbnnAFBCLdXd/QH7zSoUB7MB+14VV/FE/GxqaJfCR8F3y4wK+dM3LSfkcXVM7drra8/NIs0d7aVLhaygPTLXuNVlHmEp05jPf/2J2vlRab+4OlFshWnJfJZtium8wGJxSNXf79kq8GRhAkAJtQTpdrmt/TJBmb84mC3CY9fD3B12AtRtSIXOkd0qDtQHZD74BlRI3/TP97TbulSo9B5Hh8pOFQfojKCj3fivZVtl/sZyZR6IJqfCNSgs08rOP/dfh2xa/HaFdyzG0QFo12yX175N4v8VqMVOmb/hKJ1nvaotsAeARaWaIL1fhVFyG5DZomCbZEZjt3pfb/YdE2daIgpuVqFDxM4zd1UIwm2gYb+3KZwDKlRK7lShc6UewTXZ/RWZ/WmqQDsgUIufvcb0qvh60KNCZ6+/o88G9f7PkJ0q1NOoly1CJxVG0Pl/B6pn/256ZILvUqPjdvR8v5hOBAAVRVknfa/MTdFumVGQTplR0q0yKZ/Bdbrlvc4a2Y01ItMhsl3FN5c5b9t1Ki7StD1kX9vBUu/NKXUIsBDsk/lb4JoVv30y1x/J3MTbDrxDKgTh21X4LBlU4f/ikIozfuKYM+4vPtej4k5FV/wOAFFdp+Jpj7tVWN3FChs9HxIZdwBQUpR10m1gPiRzE2VHyiuNvtrRWy7CjTUo01li00BXq3CDW2nf4LrE/krK/mMclZ+rafcp9dha+lCgpQyqMKUH8RqUuY4EPxP2yXym+IPvkZB9bXZQHEE6HYtAPOwqQNt92/pkgvJeFaZH2o4vO5CwSSx1iOZZq/mdsRskHQhs29GEthwMvOft3sO/bbQJ7UCLq7ROur2Y+udAR3Vz5V0AoKV0ylzrdsrcaAazhPxyohJ4tcLWRC9lWJU/RwY1//8gSvG31RHbACCaQZm/2d0yAbmdiuc3pPmDA0AzjEr6tKSPldlnTM0J0rdIylTYZ1sT2oEWVylIl8yNapfKj4baNWZzojougPYUHEHvUfmaGizXBQAFQzKj6rbmjd+gikfbgWbbIRMgryzxerMC40odBk9KeqRJbUELqxSkR12Ptsvbr9rRdgBoFdWOjFNoDAAKemRG0sPWP+8X2UdI1qRMIH5fyGvNDozLdRhsaWI70MKCQbq/4m1w+6Ey5/HPNQKAdsVIDwBUz7/sq1QIyHu9h63u3iOTmcn9IpJwv0ygvj6wvdmBcakOg0+L+ejwBIP0EZl0pWDRpCgj6nb+OgAsBL0y1z07KmTX+yZjCAAMO3ruv0cclgnE7ZQgfwDfJ3NtvVnR61MAcdomUzDOSiowDnYYnFJz5sSjTYSlu9viH5IZId+paHMvh0TPKID2Z9fuLrUM15DMDSbXOwCLXTDLcrvm3y8Oylw37XW1U6ao3FZVXikIiNtBSY9KuknJB8b+DoMdMiPsgKTwID2nwrI3vTIB+rBYSg3Awtcjc/PYqcLI+YjvNZu6uV+mSCaBOgAUj56Xev06mYGffm9b2Nx1oBm2yQTpSQfGtsOgW9KuBNuBFlSpcJyt1AkAi8FOmRvHPTIjQsEg3J/aOSCm+ABA2Oh5mJxMID8kcx0FkjIq6XfUGlXUt0lalXQj0HqCQXqfClXa98mkJfUHDyqDEXcA7cqOlNv6GmGj5Hbqz16Z6yVBOoDF7DpVX6djn3dMlNWDgEZphQBdolAcSggG6XZ0aFDmImqfR2WPA4B2Y+egj6h8GrstdmTnVpLyDmCxqqWQZm8dxwLAohAM0u0o0XDgeVRccAEsdMyjBLDY2OXThlV/BtF+72vUNHkAWHSCQfo+FY+Ej4h0TgCLwz6ZUfEemVT2UllBdgoQK1oAWCzsdKBy9nv7DMkU1gQA1CiVdAMAoIXY5YB2qzBy1Ok9emXmog8E9gUAAABiU6pwXK0oHAegnW2XCcj7ZYLxUjU5totrHQAAABqgVOG4WlE4bpHYd+hPXzl58fnLLlv+5kd/77q/2Jx0e4AYbZWZ6tOrwki6vG22I5LrHAAAKPL5/+/3706n0u9q5Huk5Pxgy68+8AeNfA8kr1ThuFpROG4R+MfDf3PP+NnvXSZJJy7++KbHDu+688bubX+SdLuAGA2qcC20QTrzzwEAQEnpVPpd11z+zvWXvuayhpz/9PkTGjvxDOuqLwKVCscB85ydevF9/uevXHju/5REkI6FqmWD87Pb/kd3WnOvfljPKT15ya4/PlzLuS5su3NjYNPoil1/wvqtAABU4dLXXKbXX3Jl0s1AmwsG6UBZ+w/vvv7o2Scu8W+7OHeq45Ef/OVnNv/Sn340qXYBNeiSmXueU2HUfGcVx/uPS0RqSfoedyb/ah2R1M+tmjz/fPbpop0++9CG0IM/8oEn/U/dwH6u3Dsk7YitsQAAAIikVOE4O+/S3sRGReG4BW7WufChsO3Tc6f6JBGko53YGhx2mk+nqqvJUe/0oPp1pDq0slP69V8yz5cuXSVHxUH5bR8odXTxfr/77wvff+WJ+NoIAACAqpQqHGcLwFVbSC5K4bi13lfSKNvQmekXbwjbPjGVecMTR+5Z8+/feutYs9sE1MgG2TadvdqR8dZIg1+2RLri8vrPE8c5AAAAULdSheOGA8+jilI4bov3lTTKNnR65sWSd/IXpk99QtKtTWwOUI8RmaXU/ILPm8HOBT+YwHsDAACgxVQqHBd2E4tF6sCz997648nHnVKvX5zNBQtPAQuFXTu9U/FP66k5SHfXv+Wzztpf+HyMbTH+04anU/948P7YzwsAQAPc/e2+Lco7ayvvWbsPvnsvA4xommoLx9k569awpCG1StonGurszCs3lnt9cnpsXbPaAjTQTplr3XUy17YuSYdUWIpNkvbIrKeeqFTP2553U278J157xSSV3QGE6JTUW+a1SvsADeJscYI1WeJHkI6miRqkD3iPzpDXcjI3rIy4L3AXZid+tdzrs/mp1P7Du6/f1L31W81qExCznSrU4eiVGTHfKXPts9e6fu8xoqQLx6Etffz7maKso9lZje5817qqO0V2PJVZNTurbv+2jg4d3vGOdZP1thEooUfS/hj2AQCUESVI361Chfdh72FHzntVKC7XqRYYWULjnJ/N/XylfaZ0doskgnS0K3ut267CChd93rZB7zEic13sE0E6auDKOeB/3n3Z8rGfTZ5+NUj/3E8nSo4GffhNq19dOi9zbmbV1392Zr3/9dlZ9wZR3wAAgLZWKUjvUfFNa9gNaZ8KgTwjSwuUtz56qtJ+M7NneyrtA7SoXpnORv91zJ+yOeR93SdzzetRYYQ9Ee7nHjTBXull1mpTam11xObfXbpM/+7SpZKkSzvSayStsa99+E2ryx366v/Nv1nRoc1Xvk6SdGJqTt965XxD2grIXPeGKu5VnSjFhgFgUaoUpNsb1CGVDr79S7X1ltkPbSyfmrk5yn5nZ1++ptFtARrETufxB92202lE3FAiRq/rSOmKFUvqOseylFP3OYCI9iTdAABYTCqNjNqb1ko3p/b1sDnrWACm5s6+M8p+52ZPXNLotgBNZIN0/zWwy/uaE0UzAQAAELNKQfqI97VSCrO9aWWkaYE6P/Pym6Pu+43Dn3tfI9sCNMg+maC7RyYryNbckIrTPG1hOa53qMmGn3/NA2+5dFms57xsWYc2/PxrHujo0OFYTwwAQPNskeQGHgrZtiWR1jVRpSB9j8yNa68KN6ZB9rVhUeF9wZqcfmFV1H1ndHFTI9sCNJBN6dyvQnXiYd/2oyoUkiP9EzV526XLRl+3pGKJj6osSzl626XLRqnsDgBoY/dLerLCPmPefgtacE56cB10yYym52SWIepX8brodsQpJ3Mja5csQkzGx8e7r7zyykRHRr5x+HPvy5yt9PdSMDU7+fYGNgdoJNvR2C8zfWdIkr8eQ06FwnJc6wAAAOK1Q9KBMq8v+FF0aX6QbgvAldKlQrV3v05ve07cuMYmm81uzufzX81ms//1jW98Y2I9Rnln5qZq9p/Jn42cGg+0oO3eI6xy+1aR5g4AANAoByU9IOmWkNce1SJZZjQYpA+rvurs3LzGa5v3dbMSTOuYmjvbXc3+Z2dPrmxUW+rxhWfH7lyacg7f/NY3fjnptqAt+AN02xHZKdNZ2RKdkc6Hf/8GN+WW622uzUc+8KR2MXsJAAAkYodM/BOMKbaF7LsgBYP0fWrizWcmk1mVTqfvT6VSO5JO6W41mUxmrQrr4VY1kh23MzMvrq1m/9n8RWf/4d3Xb+re+q1GtalaDx0Z2/bcOXf7ZUvcSUkE6Shnp8zUn+tkAvUuSYdUvHrFHplRdaBqf/t87vYKa6HXfF6ZGxsAANrZqKRdkm73bbvD274oxFu5pkrpdHqjpJvy+fyuJNthHTt2bOPY2FhL3OCkUqnN/ufHjh3bmEQ7MpnMqnMzr6yo9rhz+ZOR1lVvllMz+ogknZjRqvueGf3NpNuDlrVTZspPl0yNDbvNpr4Pel/7VX5qEAAAAGq3Q6ZInCSdkgnaF41SheOGZUbUS81BL8UeF4nrut2O40jShmPHjm286qqrEp1j4LruDsdxNmQymfvXrVuXaE+N4zhFQXo+n+9WAnMwfnz6azUVZ5jNn+2tvFfzHJ9x7TKBmnWcj0r6eoLNWeh2yFxI27HKtL3ebVfhGmiruQ96jxFJu73t9UwPAgAAQGnbJH3V+9qO95U1K1U4zlYurlRILqiqiseO47yaRu267hYlWAjg2LFjG13X3SBJHR0d25T8nIcN/ieO41Q1LzwuU3NnNlfea76zM8evjrsttXr4yOjHf3BOjn1+ds79lSTbswjcLlN5c4faa4mMXpkRc1u93W6z7Frp+2SC9B6FF5drGvdzD5r56Ld9IN4Tf/ahDZV3AgAAaKhHJH1a7XU/GYtSheOGA8+jqrZwnH+uc03BYFy8TgL/94kF6V6HQXBzIkH6dP7c22o57uLc6aWPDu+49qaeHc/E3aZqTcw475UK/54nZhR5zXfUbI2k+2SC9W2S2qHmhJ1z7g+67ZKUI6IwJgAAQLMlPXCaiEqF40ZUWDe40VYmlfKeyWRWqbjM/8psNrslqWXP8vn8Rm8agN/6JNpyanr89bUe6zpLb5X00RibU5NTefetwW0PPDt66y1vX3tPEu1ZZDZIekpmKY12TFWyQbo/QLdTJ3JKcBQdAAAAC1M9heN6ZAoq7VThRrYq/nR3TyKj6el0et77vu7x73wyibZIkuM4oUXixsfHqx5N/8bhz73viSP3rKmlHU8eefA3Z/MX5/UWRDWVn/z1Wo+NPDrC0wAAIABJREFUy97nxq7PzWppcPtsXr+RRHsWsVtkKnK2RGHGEvbJBN09MmnuvSpc24Z8+9kpQIysAwAAIHZRgvROSXsl7fdt65dZkmjAe9jvq+K67prA80QqmCukcyB1+twVE3/35aSqgIfOB/WKx0X2/z6966HM2Se/9NOJfxr9h0OfeKTaRpyZful3qz3G7/zMyTfXc3wc5hz3Q2HbL7rJzPFf5FbKzFcfVcLTW8rY433dr8I1b9i3/agKheT2CKjBh67pvKOdzgsAAJormO4eZq8KxZNskST/SNKwTNC+U2a0qZ7RpURSulViHfL0hYufUZOrgJeYjy5Jcl038nrlTxy5Z01m8p/fb5+fuPjjm752eOfHf6t7e+QMgYuzubo6Tabmzizf8+2+7v5f25fYfOSTM7ohbPu06/5Cs9uCV62RqdT5pMyc9VZa89JO7+mXud4NSfIvJ5hTobBc5CKZqN7d37y5OOsi5Y5+8Nf2VT0F6e7vvH+tZvPFq1R0pO7/4Lu+1Eq/dwAALAp/eiizK9XgwbJP/vLapAZ+Y1MpSLcpnzmZG1W7PrCdk7lVhaC8X2aEKVKQ7s0Dn6fZ89Kz2WzJEb0lYy9enXvwsTWdv3/jWKl94pbP57tD5qNLqq7C+4XpU58Ipqqfm375DyVFDtLPzrxcU5q8X+fSa26VzBrlzfbwkZE1PzinN4S9RvG4lrBBUkamaucOtc589e3eI6xyu/+ahwZyHN3uf/76116pY+e+dZ99/thTnyl57I3vKJTCOHl2XN/96VeKXndndVCt1TkEAMCikHKdbrdE1jAKKgXp/vmYQ4Ft/mrHQzJBeuS1sTs6OrrDRoybvR644zglR64lSadOPSCpab0xpeajeyKPpJ+devF9wW2T0y+s2n949/Wburd+q9Lx+w/vvv7o2SeiZFqUdW72ld9TQkH6ebfjE1K+5OsPHxlZ8963djWtAwYlfUyFJdt2JdwWv7CicC0VoDsf/v0b3JR7IPYTf+QDT2pXs2qGlvamN/yK/u0vhK+Y6A/Ey3n9JVe+um9YwA4AANBqKgVhdkmiEd+2sGrHsUmlUpED0ZiUnRu7dPTFd+e2fWpV567/1pRRPsdxQjsvPJGmA3z3R49ee/jEg5eEvTals1skVQzSZ50LoXO5q3Vx7vTld3/n/WuTSC09djF/S7nXp5T+HbVWUJiEpq+mUMJKSZ9SYcm2pNvVqULHY6dMavugty0nUt0BAADQIFFHSm2w3qXKSxJF4rpuaLqx67pNK+g1Pj7enc/ny6Z0O9MzqdSaX7hfTSh0lclkVgWL6QWNj493X3nllWXneOemMiWHwKLOMz87c3xTlP2iSM25m9XkYPjvnxk7MjLllv39vjDrvrNZ7WlhrZZutF7SAUmPygTrSaQk24KYnb5tNpOoTyZw3yOT+g5U7W+fz93+4Tetbsh51dorKAAAgAgqVXe3N6a9MoF4n+81O5JkR5z8+1dUJhhvWpA+NzcXKWDtePHEb+UefKzu+dmVpNPpiu3J5/MVMw3OTL8YWixNks7OvHRVpeO/+6NHrz158fnLKu0XVUdq+cfiOlcUX35u7F9Gpty3VNov7zoV90FibpKZr97sKvx9MkUwpfDicPYaZ2twoEFuXH/bDaVS3Wv1+kuu1I3dH3vyg+/+UtKZGgAAACVFCdKHZAL0oyrcvO6RSYHv9bZ3yYyix7Ek0coYzhGJ4zhbKu9lRtO9uekNFSWLIMo+p2devLzUaxfnTi99dHjHteWOf/n8s6UrMtVgJn9h7d3ffH/D5/UP/ejYtZ99evTFZ867vxxl/7NVVMtH0z0p6R2Smr0ygD9A367iqT52+6D3fb8S5n7uwQO666H4T/zZh1otwwIAAGDRiJLufrOk3TKjRnYupk3zHFFhmaLtqiLd3XGckgFSlJTuemUymbWqYsm3ZT/Jbph46H9dv/oDv11xPnetKhSNs/uUDSwPPHvvrT+efDy8PLwnnX7dRyXdGvZaJpNZ9eRLd767UjuqtWrZFf9DUmzp5Q88O3rrkpRziSTNuc7bz866G799du7q6TI1AIOo8C7JBMNxqyfAG5NJc38kprZUo0emw9EusVbKHpl0+J4y+wAAAAA1iRKk2+XXwpYjGpF0nWorIlcy2Cw1Xz1OqVSq6jnmHaMvfU3F81TjFiW1t2yQfmE293uVTnB25qX/UOq1H576ypcuzk1WyrCoWtpZ+svZbPbw3NzcxnXr1lVdhC+TyaxNp9M7HMfZ+P0z+TXD51xJNiKvIjIP2Hsk+76b3/rGL9d8gvbXiAyHWv5DTsnULUhyPm3U2hr29U6FXxcRg8eevutA1AruVZ338KfJEgAAAC2tmmAseCPaI5Ma2qfaRpRKBuL5fL7hqdGpVGpbtcd0HJ9YNXnXFx5vRHvGx8e7FS3Vv+wN5unp8YqTOCenxq7af3j39cHt+w/vvv5n558uGcDHYH06nb6/2oOy2ezmdDp9WNItlQrrVev8XP79cZ4PNXlApvOpXQpedVXeBQAAAKhNlCC9U9JeSft92/olHVKhCrL9vholU81TqVRDR9KPHTu2sdZgb+nR8ffk/u7LH4+7Td768JF4qfqhJqdfiPRvl5sZKSqI9cSRe9a8NPXMP83mL5ZNlY/BTdlsNnIWg7fvV9WgWgWvzOo9jTiv9fCRkYYXHGxjT0q6QWbZtaYscVjBPpnOyB6V73i0BeOGxCg6AAAAYhYlSN+rwrJDNtXbBuTDKhSL26mY5mhGXYYtk8msOnbs2K5sNntwbGws8iic67qRCsaVsnRk/M9PfuGR99VzjhCRswc6OjpCg/THn/5M5OyAianMG7586I8PP3HknjUHnr331p+deeqH52ZeWRH1+Ho4jrMrk8lU7EzwsguqHnmvxuSsVtz9dLYho+l3Pz32vLcWO4qNSfqvMr/zrVZl23Ze7Vb41BbbMenfF6jKh67pvKOdzgsAAJqrUpDe6z1ykjZ5X/tVSPfc6j1soB5pSaJf/MVfrFRVO9KcwY6Ojh2u635M0gbHcW7PZrOHKwV/3ij0LVHOX4ozPeMsf+boQ3EuyxalaJxVajrA1NyZqubZ5y6Orv/pxD+N/njy8c83K0CXJNd111SabuCtGf+ImlDt/0I+f1fc57z76bHnx6fdq+M+b5s7JekOmdoLDe18qcN2mc7HHpmVK+w1rU/ShIpXuIhjNQsAAACgSKUg3Y6M26XY/NtGVCgY519PPRaVUqK9lPXg+tvr0+n0wXKBejqdjmXeqzM9k+p4fuyHcQTqXkAa+TylKrxfmD0RuVp90hzH2VYqbT+TyaxKp9MH455/XsrJWV1277OjsUxhePjIyBoC9FCPygTnO9Qaqe2l5GSKYQ5639sOyS6ZkfURmUB+a+jRTeZ8+Pdv0G0fiP/EH/lAI6r+AwAAIIJKQbpN9/SvFWyD9FoqukuSrrrqqijrU5cM0sfHx7u9UdYwJQP1Y8eObVSdo+h+6dzpFXEE6ul0utpK86H/filneTstKbYyrIicDdBVxfJ4cTh2UX9e7/zxh46MbfvxhRQBerGnZeadb5Y0mnBbotgtky20XdLVMgH7zd7jOm9bueXZAAAAgLpEre5ug/UuhQfpUZcukiQtX758eYTdbgkbac1kMqvy+fz9Kp8GvT6dTo96QbmkioF9zWygXucc9Wqr2c+bDlBNMbYWsiGbzT6SyWTWZjKZVdlsdks6nR5VkwN0SZp25YxMpX6446nKc+WDHjoytu1Th0dzz51zP3UuH2lZw8XglMy882613rzzUnplAnR/Ecxhmbnn+1RHxySqd+P6225oyHm7P0aWAAAAaGmVgnR/GnuXiuec26JJnTI3tv79y1q5cuUbouwXHGmtcpR1peu6B7LZ7OFsNnswn88/pQbNb07nTq9Y/szRL07c95U7qz029+Bja9Knzv5utcd5RdVeVc2c9hZzUzqdzqTT6Zyk+9SEOeilTM5qxWon9UyUfR86Mrbt80+PPfWXw6Mzz51zP3VipvSSgovQp2WyPVp13nkptjOSiu1oqL99Pnd7O50XAAA0V6VRPzsXvVemiJK1RyYFvlem+nunzI1tpEJKy5dHTsvekM1mD8vMY13lOM6OGuYpN2VU1pmecZY//ZPtp/7q3i1zV1zWt/oDv/2tSsfkHnxsTer4iSdSp3/+dXMrL6nq/bwl2w77NrXjSHrLeXnGveozh8eyVyyde/d739o19tCRsW1zrvP26Xz+mgt5vflCXisnZ7XiuXNu0k1tZZFXGWgxdrS8R+baFqnTMUnu5x48IEmxz0v/7EORincCAAAgflFSc2+WmafZJxOI71OhaNKITIA+JDOHM2q6ezWjpetl1smW67Z+YLTkheNvWPLC8W+e+qt7X8qvuvSx/DLnWUly5pzLnQsX3yVJbjq1vOPk6Tenn/7RSmd6xtE731rLW22UN1KZyWTWNqvI2mLw8ox71cszqdH/6/ujMsF46//eITZ7ZNLd98sE7OUC9ZyYn94wjz1914Eb3/HR+M97+NN0QAAAgJYWJUjPyQTqdrTcb0SmmBJzNQOWvHD8DXrh+P/WqPP709tTqRSj6ED9OmU6JK0eFWpwhBkWQToAAABiVk2Rq1Kj5AToCXBdd00mk1m7bt260VQqta0dsgyAFlftyDhz1wEAABC7aoL0HpmUd7sW+naZVNABmRT4kRLHoUFSqdSWbDY7Sqo7EJvtSTcAAAAAi1vUIN2uHexnKyEPeI9BkfrZVI7jUMkXaJxemc5Je60bkemYpEMSAAAADRNlnfSdMgH6iMwoUzC9fZ/MTeyACjezANCu+mRWs9gvc/2zHZG7ve12RQugJh+6pvOOdjovAABorkoj6T0yN6eSKR43rOK10iVT6b1LZtSpX4ymA2hf/SoUjxuRuebZkXNbSK5P5pq3SQnPS3c+/Ps3uCn3QOwn/sgHntSuZDP/b1x/2w2SYv/Zbuz+2JOmnwUAgPgdGf9nLUkva8i5Z+amGnLexepbr5zX0pTTkHNP5+urF1YpSLfzz+36waXs8/btFUE6gPZlOyX3qLDUpF+XTIRnOzCZww4ArW2LzLK1TXHt/X+99t9cUs1Kw9V77uTLuyRNVmzHN/6NdOlrG9eQF45L0sHGvUF0H/y1vZL0ZIPfpuLP+pf9X5997aVLn25kI86dnp6N0pZW9clfXiu1wP/VXbf8zokVqzrL/h3V68JkLsr/1f3eo0ilIN0/F7McO5pECiiAdtUrE4TnFB6gS+ZaOCgTqPeJIB0AWt0a79EURydP6ujkyUa/zfoWaYckbWjGm7SL5585kXQTENGxpw8l3QQrNIiPMiddqhx8d1XXFgBoOVE7JYe8r11KuGMyP/zcNXqlARn3oy+surDtzrXxnzi65154svvC9OlYzzkzN6WXTx9ddd+BzatiPXGVnjs9tfbMTD7Wc07lXT13emrtjqcyif5sAACgfpVG0u0Sa70V9vOnxQNAOxqSGUX//9u7++BI7vrO45+Z0dorWOwdG0ICJrsrkwceDi/R4tjFGptEW0k5dWUnlBbiJGdfQsk4wGVzXLELoeIlxdkScOUNgTKrcwWnfDhEOoLN2TkTifP6IcEPq3hN7AVsSyON5MfdlbSP0uph+v7o+Xl+05qHnpmeh+55v6q6Vprp7vn1zqjn9/09fH8m+C4W/XZn/50osU9DxJ7+yac0/ZL0e+Vu0RW676FLMu4Q0b3Bnti/yaM/vm1dR6d++Rd+PbBznlg4ooMT913idKzfqiYOFXzotTPXf+ACR5de2BnYOY+eXdFDr525PuYOmQvtMEggIPwdAAiLgvcqP0H6mNxK6X6tHQLapVwQPyeCdADhNSd3KHu/3OHsO7U2CE9mn5fIvwEArWrN/E4ACBM/66TfKHcpoj65wbi9Pnq3td+AWD8YQLgl5TZM9shdbs1eF92sYmH26VIuYDcG1ej74Nll6cVX3Z/POUd6S5Uj8M05WsjC0gkdOzUjSVqXOFfndb6l4nMsr57ViYUjkqQTZ1pnruDJlYxeXFiWJJ3XkdCb1vmdfZZzNuPo6NkVSdLRs6uBlg8AADSP35zzSbk96XaQbozKrZj67UXfe9lll918+eWX+9y9Oc59dkKJE6cb8lqL796izPkbaj7PzKkntbDa1NG3RXUmkrpowwcCOddLS45eWgrkVHV35oVnfzD24MhjnoebNoQYJSUlzdZ4jp2qfETRXkk31/i6gF/1WWsGAAAExk9PuuQO+dyZ/dkO1O01hP06cPnll1MhBdBqzHD3WpRaqrKYvaqy4ebUri9vTWj19URhq0rMb9j32UPVnGth163eJYomO/d9brKacwXh9kc+lleeRHx1vu+DwxVf27cevHbjojsH/XXrVxYP/ecP31PXZVdK+cKTqbxrW1nRZP/lWyr+v977VGrjyoryrq2jQ4f2vn9L064NAAA0TpdyPem2fkkH5Q6FBwAAAAAAddYndwioI3ceuu1g9nFHbhAPAAAAAADqxMzRdORmO/auh56UG7ibQL23oaUDgOAk5TY8eu9zhfSIhkkAAAA0gQnAD5bZb79ygTwAhFFS7n1sVqWn8Pi9LwIAAACB65dbGfUuM+TVKyqtAMLPNDiaKTz2ahZdcpejNM+Xuy+2gllF5748pNz/vXfzvldh1af8z5i5Nj+jOwAAQJswQXq5YZ19IkgHEA12Ho5xuUPbd3seC8PUHnsqUneTyxKEUkF6FL5/vMG5vc0qGu8hAH/s6aZmo7GuOXpU+run2PfRbvGetYIuue9FoTrEbPa5MNTp1jAfzFmV7qUwF+5NLAcAYeTtNTfbkMLTY3swu80qGlORzPeM9/8/qdx7FdaVRkyD+LjyryGp3HWPN6FcAJrDdH7ZG3Xs5qgmSLe3MIy6i6IelW789m6mUyZUzAWapda65FYcknJbHkwFYla0GAGIDm8laVzhuceZSoVpPS7X0BoGxYJ0KTflKoyV2C7lPl/F3iMzDSOsjRAAKlMouJhtaonaV61BOjm7Gs80fFezhaoe0S238lDqgmYV0qECAOCRVP7c9BHl7oHlksq1ChPQdincAazNT096GIeEm8agUu9PUuFpIAJQG9NwZzrI7GHv1LUbzxuk+7kX79ba2CkMdYcoKDSsfUiFv2N7lV/fC2WgbpZas9dFNzcPktoAiIpe5Qfk5kZtDzt21NrLr5m5jHbL/bjCP1y63Jz0UH2pWkyLf+iG2QGoCzufSL/yewVHmliudlVNkC6538V2oB727+Aw8Pag+x3C3q21OSBC+Z2ctDYAiAqzBJupCBX6IrYTyLVqZclU8HrLPBY25YL0sA4nNJWKMI4CABA8O7DrVfVBIoJRy/+/GclmNu7z9eN9n0pNISuk23N8aJLRdsntOfK2KvQrN1cdAMLMBOnlkrx0K5eUrRWVmp7Uqg0LfpSak+53udBWRE86AMMbaJj7nd3LF9ZRQ2FVS5BuN/6HvaG81XnzOFTTIGKPYgnF97K9HJH3xmAPf2/l4Z8A4EclN+RWrCh5W+0LbWFtyS8VpEtu40QYEyv5mZNuUMEDos0eMWSPDrKH8TJsurFqHcnAvPT68/aCVxuThqpRxV6ncUhrP5hmrnooLqZOzB9vGHtwCvHOyTCbPT837Lz5FWYVnfcP/pkVKiptKe2S+3lpxb+HERXP5G4C+LB+1qMapJvv2VJD83rFfQqIOm+AYAd03iCkHevbzRJkkB7WRvJW552LXksPeGgaUkwAXm5Yp8mMF9Y5gbWI2hq2xYL0KIyY6NLa5Ife+SfkWmgfJmgtdH8zDTnFsoG24nwlU4kr9Tca1kBW8jfcPazD+e110u3PnL3SQBSW0QNQnN3pVejv3a6/tGN9u1lqCdLtY8P63RsG3o63tuB3nl+rVlrrzfSAmHkQrdizVqlZFX4f7aX4wpq0xPwRjyi/NbNbufeQL772Ueq+Vere16r3OxPMlWqpN5XAMN6ryiWOC8X8sRIKrYtMzxnQPspNIfXOlw1rXSxsasnubr+nYfzeDQv7/YlMPT7uc79yrfft2rrfJ/faByRNKNqVqDFJw9mfw/h+98kNXkYl7ZB7PcaY9Viv+OJD+Jih+6PK/2x7DUqaU/TuVaOSdmb/DasdkvZo7fs3LGmbcvdfANHTrfwG1kJ/74Oe36N2H48SMxrPvKfDcmMFBM9bZ59oSimawLQelRtmZ3o42q2VyF57OCoZekv1pB9UeIf1m17GUu8Pc4XaS9R60gEA4eQ3MdyIz/0QHG9PeqUbuUTqq0v5/9+RiUXL9aSPZrek3BtDn9z/DLNOeq/cAL1Xbu9MO7X0mx5X07Jprj00CQdK8CYoMQFJl9zeqjDqlvsZLdXTVqoHEgAAoB7suqO3x9xm17O7RG96GHQr/B14aFH2XORi26za70ZRKItyqczKYVEucVxYW6gOqo2SScAXetIBAM3mXTqz1JQ7e9WlSM2/bWG19qTTo15fbduTLrm9i9tUeK7cmNwWv3abL2daxbw3S/NY2HvTxyTFPNsFcnuh+9V+DTIAAAD1YNepRlV6Tq131Cp5dBrvYq2tI3s3k2PEtlsRCiBbyJzn97b+e0haW7vyrsfn3cI8T6jYnHQp11oVxpZbP3PSDRoh2gM96QCAZvKujV7NRuBXX7Uswdal/NHIYR9t26rsDtOwLsW6ht/s7rY5a/Nql9aLPhXubY7Jzd7IPKHWY+ail/oyM8ubAAAA1FsQIy/DPnozyiaUn9U9CqNtW5Gdb6rWJNBJufFAJJJJd8vtpZxVe8y3MIFcsT8yk3QtrC05xXrSTfLAMLfa2uuk240o9jrp9I62D9Mjbj7z9mZavscLPGc+R3xWAAC1sNfRrmWjY6h+aulJl9aOltgfaOkg5WIze/pxtfo85wplsL5b+UtBhDl4q4S97FoxZkm6ML6x5RLHjSu8Q3W8w44KXVu7jAbB2mQ9lW4E6QCAanlX0qkk0PYmywrjNMSwqDVIl+d4gvTgef8eavk/tmPbpiac7qhw/265LQy9ygVqE3KTx40q+ktYmQQd3mQQXsPZfXsVrf+TQbnDdgpNdQiDCbkJP/rlfo7txoYBhfvaULlR1bakIJ8VAEC17BGZlS5jPCH3O8z0GJp6Od9LrcfbscV7FLwJ5WIvyf3bGlTlMViP8nvhB4rt2EoK9ZqPiHX/AAAAgEp4VweqpufPOyy3HUazNkOtPene94k56fXhHZlS6ejYUCX5s+ea24G5CdbbYf45AAAAECRv4FZNp5c30A/zykKtrNbs7vZ71NKBXwR4V98al7+/rW6tnQ7bco1eJqOdncjiYPYx86EstTQRAAAAgOLs0am1BNdmedkgEmahsGqC9C6tTWbWkoFfBHlHfptO5kL/933K5RFr+bwBJgAfz/5cKPEZQToAAABQOW+Sq1rq094AkgRywfP+H1e7hXXVpzAqFKj73VomQPeuk26SGXTJDdB7RLZrAAAAIAjeLO6VJIzzGpWbNMs+N8OpW8+gpB3NLkQb2SHpRlWWpG9CbjLhG+tSooCYuej22Pwh5RId0JMOAAAAVM6uXwexlKd3WDVDqoNVbU96v9z3gkaT5uqV+z4UWoJ5NvtcJcsftow+rV0vzsxXJ0gHAAAAACBAiTLPj0m6S9J9cocMdEl6T/a5t0nqlLQo6eV6FRAAAAAAABTn7V23s+YxpAMAAAAAgCboljvk3R7jzxB4AAAAAACazPSuE6QDAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAQCtKStotaUiS49nGs891N610AAAAAAC0gaSk/VobmBfbhiR1NaWkAAAAAABEWJfcXnK/Abrds06gDgAAAABAQJJaG6CPyB3W7tUjqd+z78HGFBMAAAAAgOjzDnHv9XFMt/ID+766lQ4AAAAAgDbRrcoDdKNH+T3vAAAAAACgBnYv+v4qjh+xji+lR2t77GflDqlPFjmmy9rXDL3vkzu83j7PfvnLNl9NGQAAAAAAaJhZ5QLWahLA9fg4rtBSbt5AuafAcXaQ3q/8BoFCW6kh99WWAQAAAACAhrCHutcr+Vu54NjevL3hXRUcawLtQj3itZQBAAAAAICG6FVtQ93L6VN+ALxf+b3uPcrvHR/3HF8oSO9XfiC+2/O8NyN9rWUAAABoiN3KH+JYy5w8O3GQfX6/iYS8ZemvogwAgMrZAWyh5dZqZc8dL5WQbqjIft4gvdhwdjtQ936H1FoGAACAuurS2rVw7W1clc9JrCVI9/aAkCUYABqnkiC93HxwbxCdtB4fKnNu79zzQo+X6uG2h+3bIwKCKAPaTLzZBQAAtB0zzG9C0g5Jsew2mH2+S/XpTSkkab3WcLYce7K/+0lEBACozZz1c9D3XHtu91iZfSesfYrNCS91jrkijwddBrQBgnQAQCP1KJe59kZJo9Zz5vc9kgY8x/Upfzj6fgWzVE1SucaBYesxya1wFat0AQCCMWH9HHRgagf9/fKfsK1YY0E13wlBlwFACHRp7R/8QVXfC1VoeKgZWlTunElPWWZVehkKAO3H3GMqSYpTrFIzrlxAXeucdCl/SOK4WAYHABrFboStNS+JXff0Jmzzu9nfUfYQ9FKJ7YrtF0QZ0GboSQ+3Hrl/wN7guVu5dRwbachTlqTcm1Sjhq0CaH2mZ2Ci5F75+5t7iBkaf4Hc4YBdCrYh0K4YdomhhgDQKMPWz0He1+3vmp3KTa8qt10csTIgZAjSwyupXPKJUbl/yOYP2gwf7VHjerLtIax7smUxN1x60wFUy+7NNqN6ZpULoIPs7Z6Qe+/akf29X9y/AKAR7CC9X8Hd2+0AuVmjo1qhDAgZgvTw6pMbqM/JbZUzN4CJ7O9jcoPlQc9x3mFG9uZsAAAgAElEQVSjQfZymxusnfzJlAkApNx8Pr9z7coNe6zHnL1R5Td2AgDqa1T5gfqI/DeSlhpVZSdi61Ppe7qdnb1cFvZKtEIZADSI6VEqNTfGy16j0d7sm0Atc9KNXuvcIyLxBYAce95goYrKfrmNiaan3J7L5/e8hp856fZxdoOAue9RSQKAxkhq7fKcJs+Sty5pVuYY0tp6rTdg984J7y9wPu9SnPb3U61z0oMoA4CQMDcxv4GzuTnY6w93K5eow9wIggjS7RvRrNygHQAM02DoTc5mV1DMPajQmrFdyt27zH0piCC9v8Bj5NQAgMbpUvFOJT9bsR51P+ure78L7DLVGqTXWgYAIVFpkL5fpW8G5jxBBOmGXeGmNx2A0aW1vSXejLZ+lqwJKrt7oZ4Y04MDAGg8b6+yn6C23PSoYvd6eys0eiqoIL2WMqDNMCc9vCqd11luv3oE0fY6x/SmAzAm5Ca53KP8NWfnso9tU34uiz0F9h2Um+AtiHXMd2bPbxvIlgMA0HgDchN57lF+fdK7zx5rv3LfBzuzW6HzDcj9TtlZTWEr0AplAFBHpoVxVoVbDs0cHvOc6Ukv16NUbU+6Oc7b88SQUQAAAADwiZ708BqW22JolmIzPeHmd7NWupmbY7JK9ijXq20neAtqPeBu6/x2YD5aYF8AAAAAACLDnoNZaPP2mtczu3uhjJwkvwAAAAAAtJUurU2qZIa6eyUL7Ovdr5bEcUnlJ6iblf81LgEAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAACNJV2Q0AAAAAAFjiTXhNgnQAKGyjpK3NLgQAAACap6PZBUC+mZmZraurq9fGYrFDv/iLv3hPs8sDoKG2StorGjIBAADaFkF6C5mZmdmayWQOxGKx87O/v/+iiy461OxyAQAAAAAaoxnD3VGE4zj3SDrf/J7JZPY2sTgAAAAAgAYjSG8R6XT6BsdxNnkeviaVSm1uSoEAAAAAAA1HkN4iYrFYwV7zeDx+baPLAgAAAABoDoL0FjAzM7O1QC+6JCkWixGkAwAAAECbIEhvAY7j3FDi6StTqdTGhhUGAAAAANA0BOmtoWRveSKRYDkmAAAAAGgDBOlNVmqou4Uh7wAAAADQBlgnvcneMPrEZ0/9xraS+8RiMXrSAZS0sOvWzZLyVoPo3Pe5A9Wc69SuL29NaPX1aTarSsxv2PfZQzUWEQAAAD4QpDdZ/MSp7eX2cRxnUyqV2rxly5bJRpQJQPhkpBtiit1sP3Z6fOqL9u+xv/n7m1WA8+nfz9/vrvtucOZPvj7CJ94Rf1zSZQEWFwAAAEUQpDdZx5H5t/vZLzsv/c46FwdA2P3eb77+ozdo16evK3hITMrf77e3S0tL7s8P/5t0+sz6YAsJAACAYgjSm2j27u9vjz9x2G9eAIJ0AOW9/a21n+MtydzP566TTtd+SgAAAPhD4rhmWsns9Lsr89IBAAAAIProSW+i+PFTV/jdt1Xmpc/MzGy96KKLSCAFtJi4dMC55JdvkFRutYjKvKtLemHqgUDPCQAoqFAS0KCFORno559IbY3HtLH8ntXJOJq/5dItofy/QbQQpDdR/PTCOyvZv9nz0h985lv9D6Q/v/vhZ/7+v37ovb9/W7PKAWCtzn2fO3DmhfSk6hCkx37l4ge0K9CzAgAKyEh7Y4pdX8/XiCvztKStfva9/ZGde2MxFUw6GpRPbB+K+d03Hovtc6Qr61WWeEwPyZ1iCjRV5Ie7p9Ppa6emplqyetnxyrENFR7S1JvGscUXPrGSWdTLC/92SzPLAQAAAABRFekgPZVKbZb0vVgsdtv09HRLtYodG37g45Ue08x56T88fMemo4vPny9Jc2en1j92+B8/1qyy2L57eGLTdw9PBNtzWKX/9czUrX/x5KRz9+HWbBQCAABA3fVJcmrcuhpearSUSAfp8Xj8WvOz4zgtFTjFZucrDtLNvPR6lKeclZXFz9i/v7J4+C+bUQ6vibPxn7y2nPhhs8uRSqU2Ti85n5WkhVXn2nL7AwAAAEAhkQ7SY7HYDdav16RSqbolmqhUx9zJ91RzXHZeesOdWTqWl+Ru/mz6V5pRDttdz6bvmF9R58ySc/HtT6eb2rP/6JnEV0+tun9PJ1Z1STPLAgAAACC8Ips4LhuQ5wVLiUTiWrXAWuOvjj78vo77/rXS+ehGU5LHnV09mZfkbnF1Pv6DQ9/46G9t/eQ/NLosxvRixkqs4nxJ0neaVZZXljLXNeu1gdd9/e4r9emAP4p/c7cc6UFJvhP7BO32Rz52lZTJa6C86YqhvdWda2f+cXFn8qYPDjf9ewEAEBmD2a2QPkn7sz/vkDTakBIhdCIbpBfpcW6JIH3dS7MVD3U3mjUvfX45vaZRwYmt/KmkpgTp//Ds1Nd+fMZ5/fM7s+RcPPTs1Pad79n0aKPLMvTs1Panzzid5vejy/VbGgRoT5mrvNmFp08/mvf7/U99reCRv/P+/6JS+zmKPaQW+F4AAAAwIhukx2KxqxzH8T58TSqV2rhly5b5ZpTJSBw7fnW1x9ayXvr09PRVLxwfeeeH3/vHd1Ry3IPP/O3Hfza/dpnk+cXpX6u0DEF5adn5I+9jxzP6gqTfbnRZzmScz5TfC0CtvAG33+eK7ffY89/V0dMv1lwuAACAIEV2TrrjOAXXf8wOeW+qxOzxLbUcbyfE8yOVSm28/9BXj9w/9ZkHfzb/wP984Om/fqiS4xdXT/xGocfnl9MbmjHP/+7DU7sK9VYfWXJ+s9FlkaTjq+rxPrb3qdbJfwAAAIDQ2S830/t+udneR5TL/j4rqduz35reSQ+zz/4S+3RJ2q212eb7rddDA0Q2SJd0ZZHHbyjyeEO8ds8/Xx0/eaam//dKh7z/7MT3750+9cSbze+TJ//lQw8d/rbv3vyzq6cKNnhI0gsnRz5VSVmCMLesvyj0+OmMOoYOpz/ayLKkUqmNry1rzVSAXz433tTPGQAAACLBBOh2p1BS0lzAr9MnaVxuQO61W9LB7L9ogEgG6TMzM0WDSklXNmsZM0nqmD/1hwGc5hq/O6ZSqY0zp8eu8D4+t5j6qt9zLK4c+4Viz51aevkjfs8ThOnp6ascOW8u9vzxZWdPI8vzowWCcaDebrpiaK/f4eyVuOyXPqKbPjjUlDwfAAD41CM3UN8jN4nrBdmfJwJ8jV7lethH5Sa1i2W3bZKGs8/1i0C9ISIZpGcymZJBeCKR2Neosqx57WMndgRxnnQ67WvI+3Mn7/+rlczimqzMRxefe5ffoerzSy8W3e/U6msNXYrNcZySGZ2PrDjvbVRZpOJroi+sqquR5QAAAEBkDUoayP48Z/0chKRyAfqw1madH5O003rN3dljUEeRDNKLzUe3XJNOpxveAzr37fs3rZt+pWgvcIV8lX9hZX7NfGlJWsks6rmT9/9VuePLDYs/vXyk84eH79jkpyy1yr5nxaYxuOXJqOPvnpmsOnt+pVacWJFg3Hlfo8oAAACASKvnUm19yg2fv7HEfgPZfZLZY1BHkQzSfc7Z/tb09HRDe9Sd5bO3Bni6a/wM258/my7a031q+WjZ3vjFlbmy2dIXl059otw+tcr2+vt6v86sxho25H0x49CSiNbxqesqSgrpy6evU+yTf/DhwM8LAAD8qmeQbhLCjan0PPc5qxwkkauzSAbpksr1pEuSHMf5s+np6clG9KqnUqmN54zPBJrULJFIlBz6/dDhb1+9uDpf9D2eOzvxjnJD3hdXjpftEV5aOV73Zc+yUxTO97Pva8tOV6OyqxdKGidJJ1Z1SSNeH2gXxdZBr/Wc33x0Z7lsuAAANFvQSeJsZlRoj9ZmdfduvZ5jUCeRWyd9ZmZmayaT8RXMSe6645K+lU6nvyXJ9EJNdj5xeL2z/pzF1Y1veGbx7W/939WsS24c+f7In6+79+Ev1ZrVvYDrU6nU3mJlW1g+WjKxgzXkvWhGpqXM6feUK0S956VnG1Gu97v/kqPYr56T+Kqkug57H3p2avvTZ6jfAwAAILSqGRXKSNI6i1yQvrq6elUstiZPml9mvvOV50y+pHUzr0mS3nTeG79y9O/vGzl52Xt2btmyZb6SE87e+Y93r/9/Y79fbYHKSSQS96jAyIFUKrXxoVdu3V7u+GJz1o3Ty0cvKHeO08tHOh/76b3vu+xXr/lxuX0rlW10qXhawitnMx9VnYP0lZi2FXuu0DruAAAAQIsalpsgDi0gcsPd4/H4rsDPeeK03vD4Mzsu+N5Dk6+OPuw7IdiR74/8+bk/SdUtQM+6JJ1O3+l98Gcnvn9vqaHuRqk566lUaqOfc0jS8aXpwAPibIB+QD6HudvmV7Xh7sNTgX8WbMUyuwMAAAAtpFTP91j2X4awt5CG96RfffXVV506dWr+4YcfDvzc09PTV2WHr9fFuvQr58fOLP7bS48+/uG3bf/1R0rtm0qlNp777R98JXZ2uV7FsV2fTqc3x+PxXcvLy5OTZw58c/z4gx/yc+Di6nz8B4e+8dHf2vrJf/A+97MT/8f3XP1yPfKVqiVAN15dcm6Wz2Rzd/775PfmVnXV+ljs2M+tW/3Nj7y7a6rcMQsZlRzm/3fPTH78+vduvsNncQG0oS88OXnAsVatWJ+IP/3xLcm8Bsavv3DswULHfuqdF3641H4xOV/80ge2lMxdAgCIFJOl3atUPX1M7lzz7ux+pZLUjcsN5gdVOhM8atSU4e7d3d3X3HXXXXs3bdoUWOUhmwBtTY9y0DqOzife9MOxH8x9+/53Jf/gd4oGcuc/+NRjHUfnE/Uuj+XKB9Kff2pdvFMLK5XllnBiK38qaU2QnnGWfSeEO7X86sUVvWgJQQTokjvk/O7DU7uue/emooH60LNT26eWnH9+flGd7iPOxkUnfkg+5tocW9ZbSz2/Kv18hUUGavP1u6/Up68L9px/c7cc6UFJVc8jQmlv6ojrV887V5J03rr4JU7MyQu2P/lLhWcdOcrf7wMXdL7+85OzC4GXEwDQkiasn3vlBtBepfJUDUvqz/7cL3ed9ELBxG7letuHKywjKtS04e6xWOzmdDp9wM8yYn4kEol99exFz3utuROdHRPp51569PErCj1/5J8O/OX6Hz9X12RqhaxkFisO0CVpfnH61wo9vrB68t1+z3Fq+dVzyq2p7kdQAbrx4lnnlu8enij4uRg6nP7oTxach+dX1Gk/fnRZG29/euqFUuf97uHJLyw5pYMWxynZagmgAr/z/qL5LWs65ye2DzW98eFN6+K69MJOXXph5+vBejXMOS69sLP8zgCAqLAD5n7lMrBLbs/4iNxe8mJBwoQks3xxt6SDyl8HPSlpv3KB/IDquyQc1Pw56VcmEolUOp2+s5ZgfWpqapcqyP4dhMSx4+dsOPDUg6/dO/p1+/Fjww98/A0HxkI1vHB+Ob3hsZ/eu2au/cLKsZI9xV5nlo/8YS3lSKfTNwQZoEvS/Io6X1tO/ND7+J3/Pvm9p09nvlMs0J5Zci4uFainl/SZcq+9lIldVFlpAQAAgIrYQXZS0pByS6aNyA3Ub1TpZdwGrHN0yQ3KzTlmlQvah639UEetkt39+kQicX06nX7acZw7M5nMPX6XPJuamtoVi8Vuq3cBC+k4Op/Y8ODBT54YuONPVi4872exldU3nXPw2a7YUkPmoQdq7mxqj6TXx8k+9tN733fo6LfPqeQcJ5Ze/N1qXjuVSm1MJBJ3SrqmmuPLmVlyLv7KU5Nnfv6c+N2SdHQ585HnF8tnXzeB+k2XbHqn/fhdz6bv+OmZTNnjFx3nwupLDaAdZGLOobec23Fl+T0r87bODr14ZrnqpUMBAKEyIHdueZ/ye9IHs9uYSg95N+cYlRvU93ueG8w+xzD3BmmVIN24JBaL3ZZIJG5Lp9PHJR1yHOeAJMXj8QOS9I53vOOA2Tmb1byhPeiFdLx8dH3Hy0cvaXY5anF6+UjeUPVqsrXPnZ1a/0+H/scDV2/9jO+57NlGlr0KsPe8kPkVdc6vZP6k0uNmlpyLbzs0OfeWjtg3z4np/pMZ50s/PZPxVaFmGTYgOPc/9bXAh7zf/9TX9M1HdzrNHPJ+y7Ytu16aP/FnQZ/3dy86T2/beF7d87QAACpigmY/blRlydlGVXoYup/8UWPZbaCC10UdtFqQbjtf0pWxWOxKSXIc52ZJSqfTknQ8FovNN2oOejs4uvj8+fZa59Vma0+fevy37vjXP1rtTCSPbFj31p9K0hvXvWU8Eev40TvP3/GCJDmOs1HStdmtrsF5EI4ua+PRZWePqhjeM/Ts1Pad79n0aB2KBQAAIiam2HwDXmbR/67xA46TqV9JKpRxnF3xWP06QTKOGvH/D5TVykF6Kec7jlPX4G7dzGv1PH1LOr700l/JDZxryta+kjkbP5l55a0nl18xc9qvvODci//44vPaL4/aSkzbJBGkAyjoL56cdL7x/GzRDO7V+sbzs4GeD0BjvGHfnl2SdpXdsUFuuuI7ByQdKLtjg9xy6ZZDzS4D0AhhDdJRBy+ffuo/plKpja8uPPXbh2a/U9F8dAAAAABA7QjS8brF1fn4s8e/e2BxZS6QZfEAAAAAAJVp9hJsaDEzp5645Oji8y0/TxxAAZ+67qHAz/np6xT75B98OPDzAgAAoCB60gEALS/ozO7mnO944/amZXYHAAAohJ50AACa5L9/YHMs6KRxkgJPRAcAABqHIB0AAAAAgBZBkA4AEbCw69bNmnwx+LVjj8wpM/bsOwM/bwVu/9HHNh87NRP4eU8sHNHgv/RuDfzEFfj8E6mtR8+uBn7eFxeWtedHKZKAAgAQQgTpABABGekG3ffQJYGf+JExxZ7+yacCP28lVjI3PPb8PwZ+2sMzD2tVsX2Bn7gC8Vhs3yNHTgd+3ntmTqqjQzcEfmIAAFB3BOkAAAAAALQIsrsDQJS8+Gru57e/tbpzHJmTlpbcn88u116mgNhD3i/ccFHN51hePVtzmYKwlHH04oL7/3xuPK43n5uo6jzmHAAAABW5+uqrD0hy2NjaeLtKQGFXqfmfT7b22bgXAQDQglgfFgAiYGHXrZsl5SUK69z3uQPVnOvUri9vTWj19SR0q0rMb9j32UM1FrFqt//oY5u1kn9tN13xnaqu7fZHPpYXmCbiq/N9Hxxu2rV9/onU1nhMr/9fZxzN33LplqrK84UnU3nXtrKiyf7Lt0zWWkYAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAALVI3cpm/5mFyQgsyq8VM+spN1NLFeQdks6qPxri8r7BxQzIveznmx2QQKwW8WXFRuR1N28ogWmW9J+rb223mYWCgAAIAyG5FaexptdkIAUC9LNtr95RatZl/KDc+92UNEIYACvbuU+51FobCsVpJuGt66mla525a4vCu8hAAAognXSa5OUG5yPye1R3yNpoKklqt2spAlJ2zyPd8ttkOiSdHF2n7A5KPc6RuW+V2PZx7vl9qT3SBqWtLMppQPqp19Sn9y/26Tcv+Ew2y33mnbK/Zu17Zd7rYOSbmxwuYLQK/deK639TjHXLbn36DEBiJSFXbdulrS5nq+xqsT8hn2fPeRn38F/6d26mklsrGd51KHJmy7/zqSfXb84NnUwHot11KsoGcdZObua+fgtl27x9f8D1EvdPuRtok9uhXdAbvDaq/AH6cWMya0M71Y4e5v7lAvQd3ieG8s+dlDue9ilcDZCAMX0yf2cj8oN8nq1NriNij3K3ZvDyAThO+S+X7YBSXPKNUSEsRECQAmZn3vz3tjC4vU67431eYGzy4qfPvO0pK1+du88J3lPIt6xaV3i3LoU58zSCZ05e/KLkvb62X8p43Rfe9GGupRFkg68eubUciZT30YJwAeC9NqYnqlR5QLYHq2tWEVBt9xrm1A4e2/MHNVSjSg3KpzXBpRiGtaGlbtP9Sm6QboZCh7Gv+Wu7Daq4t8jg8p97wCIond1Sb/+H+pz7hdflf7voxUd8u6LPqQLN1xUl+I89/Ljeu7lxys65u2d6+pSFklaF4+t1O3kQAXizS5AiJke18Hs76bC29ec4gTKnr9qz9fuUniHgnfL7YEqVbENY6UeKKdXblA3qNzfQI+ikVzN5ATxztceVjhHNZn3pFwAToAOAECEEaRXr09uhdcE6WYoaY/CO8yynKTc6wMQDiYYt3vNzT0rylnCkwp34jgAQPh5V+gotY0oXElBx1U4oXSXop3ktNh1B44gvTpm6HdS+dnQzWNh700fk5tU0N4uUP58VgCtz9yL7GzhI57nwmyn1t6rblQu0SUAAGHQI7eOPatojHSrhz5FZzWtsgjSq1MuSI1C5ddrTrkkRWEM0sfkfyRAGK8P8Cr3eU8qmq3cg3JHDnQrfBUdM+WG+xQAtKek3EbmqI7KrdZ+hXsZ6IoRpFfHZEr29uDElJ/pHa3DzOEsFZSY3kYgCkyG80K9zRfIbXjjPtVaTGLOHhUP1Lvk5giJYgMLAETRhArHDPa2x9q/S+Ht8LOvNYy5YVoGQXrlTKbkwSLPRymBnC2pXAtWGBOsDStX+R1RfnDSnX2sX7ml5oCwM6tPFPo8mwRyZupOlPTJ/fueUzjvVaZSY+5Jtt1yA/RuUfkBgCgZUH5yZhrR2xxBeuVMxbdYkG6CvDBnTy6U3X1WuSXYil17q9spt/w9ys8KfVC5awtr9nrA1qf81ScKGbD2DatC2d33y21UDGsQa2emt3MJOHKDdnNtNCYCQLSYDiUpvDEEAkKQXhnvsmvFmMpT1FrBBiXtkNtDFUYTki6WW8H1XsOApG3ZfYCwM6tPlApUx5QbXRKluW9jcocNhjVIl9zy79TapdZG5eYG2bPmCABAFBSrY5tM8fvlxiIjyu9I8wb1XVrb0Gsae/00AJi8NXaC7IMq37DvN7v7bs81OHIb3r2xk7kG87rlzt+s6w5cR6NfMOSG5c6xCGq/VnRBswvQAHtEJRfRti3g/VrNgMIdhPsxLHrLAaDdlFs+1ATo9n5J5Qf3fSqeZG13divVmG1GnHob8Luz562ll79Q+Y3e7DaoXLLqSrTydVeMnnQAAAAAaC4zYlcqnlOlJ7vPHuUSwe5RbiRor3KB6qjcEbAmkds25Rp/+1W8J9oEqmYaqDneTBs1U+oqlVR+gD6QLX9M7khXM3qszyrbQPZ5M4q5WGK6Vr7uqhCkAwAAAEBzdMsNHoesx0pNrR1ULkC1p7bZSZ6H5Qaq9rSpMbkBp533xNtrbB6byx5vj+gaVm1TQ+0g14xqNSMATGBszl2obMW0+nVXhSAdAAAAAOrDnkddaPMuqzmm0kG6N1+JYZZenVPp4eImN1NS+XOt7d8HVTgoLZfvphRzjcXOMafcdSflf/WZVr/uqhCkAwAAAEDzmV7fUooF6WbO9JhKJ3k2y7Dax0j5QXGx15Cqy5fSpVzvdalzm+HtsQpep5Wvu2okjgMAAACA5hmQ/2WOiwWiZih5j9weej+6ivxcamj3nHI90n7Z5y42375arXzdVSNIBwAAAID6MEsA11s1wWO1AWfDglUfInndDHcHAAAAgGgwS0H72RrReNAokbpugnQAAAAACDczjLzaZcLsod6l1gRPVvEa9hD9cuc2CfWGSuxna+XrrhpBerT1yV2PsFQGSQBoJu5TAADUzgSr3SqfGX1c7vftfusxO2laqeP9Zl232Und/J7b79z1Vr7uqhGkR9eI3A+g9wNl1mIcaXiJACAf9ykAAIJhZx/vV/G507uV6xG2j7GXQOtT4V7lpKpvRLfXcy92DntpNL/Z1Fv9uqtCkB5NptI7IWmP8udg3Ci3xalH9FQBaB7uUwAABMd8n0puoHlQa9cD3y83kJXcoNm75JhZS1xaO6qtR27jeakh4aXYa5D3W+WQ3OB5SLlG+z3KH4Y+Z+1nXt/82+rXXRWC9OjplvvBHJO0TblWK2NQ0g65H+g+z3O9yh92WugDOavch3fW2rdQRXq3csNKzNySVskECaB5uE8BABC8AeUC1i65wan5fptV7jt12NrPNqHc96/kBrbe79sBlV6qrJi57LnN8PTd1rnH5X6/S24dwFsvsIPqg9lj7PpBK193VViCLXpMC5TdIuRl/kjsD1qv1iZoMC1H27Q2qYK3Umy3TpnfvRXiXrl/ONuKFx9AG+A+BQCoyhPHFup27tMrmfV1O3njmJ7iHuX3VktuADyq0kPJx+RmP9+d3ZLW44PZrbfwoWVNyP1+3Z0tnz3dbdjavEbljrKzh7N7v+Nb+borFmvUC6FhRuR+OC9Q8cpvIeNyP4wDylVg++R+yMfkVpYltzUqKfeDbh4zFV3zmBlqMpY917Bnvz1a20IGoH1wnwKAAs584a/36dTCn9X1RTo6Hn/DV//bZX52/dvH/tOhpZXFS+pZHMfRF2+6Ymivn32/8uP0NzsUr1sgvaLM4rHTq/39l2+ZrNdrAGhPB+VWUL3sISX2ZjIhFhsKaoZ7mEQLs0XOb4aX2q9VKAviQbkVbQDti/sUAABAEQx3j55KeqUMU7H1JnHw7mOGkpabj2HOVywzczVlBBAd3KcAAACKIHFc9IzKHebpnTMxoPzsyfYwTj9JkipJpFRu32SF5wMQLdynAAAAiiBIjx6T/dBOeFCOybLoXQbJ3vyuVWifb0eJ89FLBbQv7lMAAABoK2bZgXEVzlxsL0tgMiOOK3+JAimXWMkkYZLy53Ta7Me7s+c+qPzMi30ljgfQXrhPAQAAoK3Y6wgX2+yKbl+J/ewKtJ/Kr5RfwfZuhRI1AWg/3KcAAADQVsyawnbFc1zFky71Kr/SXKiHy2/lV9ljx63zjaiB6wsCCAXuUwAAAAAAAAAAAAAAAE/IIcwAAAAWSURBVAAAAAAAAAAAAAAAAAAA+PD/AeNiW6SstofpAAAAAElFTkSuQmCC"}}},{"cell_type":"markdown","source":"## 1.2 Evaluation（評価）\n提出した成果物の順位付けには、ピアソン相関係数を使用します。Multiomeデータセットの各観察について、正解の分かった遺伝子発現と予測される遺伝子発現の相関を計算します。CITEseqデータセットの各観察について、正解の分かった表面たんぱく質レベルと予測される表面たんぱく質レベルの間の相関を計算します。総合点数は各サンプルの相関スコアの平均値です。あるサンプルの予測値がすべて同じであった場合、そのサンプルの相関は-1.0として点数化されます。\n\n> **_NOTE:_**  ピアソン相関係数$$\\rho=\\frac{\\sum_i (x_i-{\\bar x})(y_i-{\\bar y})}{\\sqrt{\\sum_i(x_i-{\\bar x})^2 \\sum_j(y_j-{\\bar y})^2 }}$$\n\n\n\n### Submission File（提出ファイル）\n評価セット内の各`id`について対応する`row_id`の`target`を予測する必要があります。提出物はヘッダーを含み、以下のフォーマットである必要があります。\n\n```\nrow_id,target\n0,0.0\n1,0.0\n2,0.0\n3,0.0\n...\n```\n\n提出ファイルはテストセットの観察のサブセットのみを含む必要があります。含まれるべき特定のID については、データ説明を参照してください。","metadata":{}},{"cell_type":"markdown","source":"## 1.3 Timeline\n日本時間  \n開始日　　　　　　　：2022年 8月16日  \nエントリー締め切り　：2022年11月9日 8:59  \nチームマージ締め切り：2022年11月9日 8:59  \n提出締め切り　　　　：2022年11月16日 8:59","metadata":{}},{"cell_type":"markdown","source":"## 1.4 Prizes\nこのコンペティションは、NeurIPS 2022 コンペティショントラックの一部です。優秀なメソッドの賞金は次のとおりです。\n\n* 1位 - 15,000ドル\n* 2位 - 7,000ドル\n* 3位 - 3,000ドル\n\n賞金に加えて、受賞者は 11 月にバーチャルで開催される NeurIPS 2022 コンペティション ワークショップに招待されます。","metadata":{}},{"cell_type":"markdown","source":"# 2. Data","metadata":{}},{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nimport os\nimport scipy\nimport time\n\n! pip install -q tables  # needed for loading HDF files","metadata":{"execution":{"iopub.status.busy":"2022-09-15T08:37:18.780547Z","iopub.execute_input":"2022-09-15T08:37:18.781247Z","iopub.status.idle":"2022-09-15T08:37:34.781553Z","shell.execute_reply.started":"2022-09-15T08:37:18.781127Z","shell.execute_reply":"2022-09-15T08:37:34.779728Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Data Description\nこのコンペのデータセットは、4人の健康なヒトドナーから分離した動員末梢CD34+造血幹細胞および前駆細胞（HSPCs）から収集した単一細胞マルチオミクスデータから構成されています。細胞についての詳細は、ベンダーのウェブサイトに掲載されています。\n\n\n測定は10日間にわたり5つの時点で行われ、各時点の各プレートから2つの単一細胞アッセイ手法で測定するために細胞を採取ました。1つ目は10x Chromium Single Cell Multiome ATAC + Gene Expression technology（Multiome）、2つ目は10x Genomics Single Cell Gene Expression with Feature Barcoding technologyで、TotalSeq™-B Human Universal Cocktail, V1.0 (CITEseq) を使用したものです。\n\n\nこのデータ型を扱ったことがない方のために、この説明の下にいくつかのリンクを掲載しました。\n\n各アッセイ技術では、2つのモダリティを測定します。**Multiomeキットはクロマチンアクセシビリティ（DNA）と遺伝子発現（RNA）を測定し、CITEseqキットは遺伝子発現（RNA）と表面タンパク質のレベルを測定**します。\n\n> **_NOTE:_** 私の理解では、クロマチンアクセシビリティとはＤＮＡの読み取り可能性を表す。クロマチンはDNAとたんぱく質で構成されており、たんぱく質が結合することによってDNAは折りたたまれ読み取りが制御されることで不必要な配列を転写しないようにしている。\n\n> **_NOTE:_** 遺伝子発現（RNA）はRNAの量の指標を指していると思われる。RNAはDNAからの転写で生成されるのでクロマチンアクセシビリティと関係があると推測される。\n\n> **_NOTE:_** 表面たんぱく質レベルはおそらく膜タンパク質の量の指標。どのたんぱく質を測定しているのかなどはまだ理解していない。たんぱく質はRNAから合成されるので、遺伝子発現と表面たんぱく質は関係があると考えられる。\n\n分子生物学のセントラルドグマに従うと DNA→RNA→タンパク質という分子生物学のセントラルドグマに従って、以下のようなタスクがあります。\n\n* Multiomeサンプルについては、DNAからRNAを予測\n* CITEseqサンプルについては、RNAからたんぱく質量を予測\n\n<br>\n<br>\n<br>","metadata":{}},{"cell_type":"markdown","source":"### Cell Types 細胞の種類\n解析の参考のため、以下の論文の情報を用いて、RNA 遺伝子発現に基づく予備的な細胞タイプのアノテーションを行った（https://www.nature.com/articles/ncb3493）。\n細胞タイプのアノテーションは不正確な技術であり、連続データに離散的なラベルを割り当てることは限界があることに注意してください。これらのラベルを予測に使用する必要はなく、主に探索的な解析のガイドとして提供されています。このデータでは、以下のような細胞タイプがある。\n\n* HSC = Hematoploetic Stem Cell（造血幹細胞）\n* EryP = Erythrocyte Progenitor（赤血球前駆体）\n* NeuP = Neutrophil Progenitor（好中球前駆体）\n* MasP = Mast Cell Progenitor（マスト細胞前駆体）\n* MkP = Megakaryocyte Progenitor（巨核球前駆体）\n* BP = B-Cell Progenitor（B細胞前駆体）\n* MoP = Monocyte Progenitor（単球前駆体）\n* hidden\n\n> **_NOTE:_** ○○前駆体（前駆細胞）は造血幹細胞から○○細胞に分化する途中の状態の細胞です。細胞の時間変化の不連続な状態に名前を付けているだけなので予測に使う必要はないとされています。しかし、一時的な分類であっても問題を解くにあたって情報を整理することは役に立つかもしれないので、造血幹細胞の分化の略式系統樹を記しておきます。  \n<span style=\"color: blue; \">HSC</span>\n――> 骨髄系前駆細胞 -> <span style=\"color: blue; \">EryP</span>  \n　　　|　　　　　　　　　-> <span style=\"color: blue; \">MasP</span>  \n　　　|　　　　　　　　　-> <span style=\"color: blue; \">MkP</span>  \n　　　|　　　　　　　　　-> 骨髄芽球 ――> <span style=\"color: blue; \">NeuP</span>  \n　　　|　　　　　　　　　　　　　　　――> <span style=\"color: blue; \">MoP</span>  \n 　　　―> リンパ系前駆細胞 -> <span style=\"color: blue; \">BP</span> \n    \n<br>\n<br>\n<br>\n","metadata":{}},{"cell_type":"markdown","source":"## File and Field Descriptions（ファイルおよびフィールドの説明）\n**metadata.csv**\n* `cell_id` - 観測された各細胞の識別子。\n* `donor` - 4つの細胞ドナーの識別子。\n* `day` - 観察が行われた実験の日。\n* `technology` - citeseqまたはmultiomeのどちらか。\n* `cell_type` - 上記細胞タイプの1つ、または非表示。\n","metadata":{}},{"cell_type":"code","source":"metadata = pd.read_csv('../input/open-problems-multimodal/metadata.csv').set_index('cell_id')\nprint('length of metadata :',len(metadata))\nprint('# of unique cell_id:',len(metadata.index.unique()))\nmetadata.head()","metadata":{"execution":{"iopub.status.busy":"2022-09-15T08:45:29.173587Z","iopub.execute_input":"2022-09-15T08:45:29.174019Z","iopub.status.idle":"2022-09-15T08:45:29.549232Z","shell.execute_reply.started":"2022-09-15T08:45:29.173976Z","shell.execute_reply":"2022-09-15T08:45:29.547752Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<br>\n<br>\n<br>\n\n\n実験結果はいくつかの大きなアレイに格納されており、HDF5フォーマットで提供されています。\n\n\n### Multiome\n* **train/test_multi_inputs.h5**  \nATAC-seqピーク数をデフォルトのlog(TF) * log(IDF) 出力（DNA）を用いてTF-IDF変換したもの。行は細胞、列はアクセスレベル(DNA)を測定したゲノム上の位置（ここでは10x References - 2020-A (July 7, 2020) で提供した参照ゲノムGRCh38のゲノム座標で特定）に相当する。\n\n> **_NOTE:_** TF-IDF変換：統計量への変換。この文脈ではDNA配列において特定の配列がいかに重要なのかを示す指標への変換。（ただしDNA配列の区切りがどこにあるのか理解していません）TFは特定配列の出現頻度、IDFは各区切りで特定配列が出現する割合で$log(TF) * log(IDF)$で表す。\n\n> **_NOTE:_** ゲノム座標：ゲノム座標は参照ゲノムの染色体名や開始位置、および終了位置が次の形式で含まれたもの。`chr1:1234570-1234870`\n\n","metadata":{}},{"cell_type":"code","source":"df_multi_in = pd.read_hdf('../input/open-problems-multimodal/train_multi_inputs.h5', start=0,stop=30).astype(np.float16)\ndf_multi_in.head()","metadata":{"execution":{"iopub.status.busy":"2022-09-01T01:26:23.836854Z","iopub.execute_input":"2022-09-01T01:26:23.837238Z","iopub.status.idle":"2022-09-01T01:26:24.915313Z","shell.execute_reply.started":"2022-09-01T01:26:23.837207Z","shell.execute_reply":"2022-09-01T01:26:24.913871Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"\n* **train_multi_targets.h5**  \nRNA遺伝子発現レベル（同一細胞のライブラリサイズ正規化およびlog1p変換したカウント）。\n\n> **_NOTE:_** ライブラリサイズ正規化：同一細胞内の総和が１になるような正規化。\n\n> **_NOTE:_** log1p変換：$f(x) = log(x+1)$","metadata":{}},{"cell_type":"code","source":"df_multi_tg = pd.read_hdf('../input/open-problems-multimodal/train_multi_targets.h5', start=0,stop=30)\ndf_multi_tg.head()","metadata":{"execution":{"iopub.status.busy":"2022-08-30T05:04:38.344971Z","iopub.execute_input":"2022-08-30T05:04:38.345457Z","iopub.status.idle":"2022-08-30T05:04:38.496414Z","shell.execute_reply.started":"2022-08-30T05:04:38.345417Z","shell.execute_reply":"2022-08-30T05:04:38.494878Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<br>\n<br>\n<br>\n\n\n### CITEseq\n* **train/test_cite_inputs.h5**  \nRNAのライブラリサイズ正規化およびlog1p変換されたカウント（遺伝子発現量）で、行は細胞、列は`{gene_name}_{gene_ensemble-ids}`で与えられた遺伝子に対応する。\n\n> **_NOTE:_** gene_ensemble-idsがよくわからない。gene_ensembleとは何か、なぜgene_ensembleで区別しているのか。\n","metadata":{}},{"cell_type":"code","source":"df = pd.read_hdf('../input/open-problems-multimodal/train_cite_inputs.h5', start=0,stop=30)\ndf.head()","metadata":{"execution":{"iopub.status.busy":"2022-08-30T05:04:43.768712Z","iopub.execute_input":"2022-08-30T05:04:43.769149Z","iopub.status.idle":"2022-08-30T05:04:43.947673Z","shell.execute_reply.started":"2022-08-30T05:04:43.769115Z","shell.execute_reply":"2022-08-30T05:04:43.946380Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"\n* **train_cite_targets.h5**  \ndsbで正規化された同じ細胞の表面タンパク質レベル。\n\n> **_NOTE:_** dsb：たんぱく質信号の背景処理","metadata":{}},{"cell_type":"code","source":"df = pd.read_hdf('../input/open-problems-multimodal/train_cite_targets.h5', start=0,stop=30)\ndf.head()","metadata":{"execution":{"iopub.status.busy":"2022-08-30T05:04:53.034062Z","iopub.execute_input":"2022-08-30T05:04:53.034555Z","iopub.status.idle":"2022-08-30T05:04:53.093595Z","shell.execute_reply.started":"2022-08-30T05:04:53.034514Z","shell.execute_reply":"2022-08-30T05:04:53.092348Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<br>\n<br>\n<br>\n\n## Splits　（データ分割）\nデータの分割は以下のように行われる。\n\n* トレーニングセットはIDが13176, 31800, 32606のドナーのみから構成されます。公開テストセットは、ドナーIDが27678からのサンプルだけからなります。プライベートテストセットは、4人全てのドナーからのサンプルで構成されています。\n\n* Multiomeサンプルについては、トレーニングセットは2日目、3日目、4日目、7日目のサンプルのみから構成されています。公開テストセットは2日目、3日目、7日目のみのサンプルから構成されています。プライベートテストセットは、10日目のデータのみから構成されています。\n\n* CITEseqサンプルについては、トレーニングセットは2日目、3日目、4日目のサンプルのみから構成されています。公開テストセットも2日目、3日目、4日目のみのサンプルで構成される。プライベートテストセットは、7日目のサンプルだけから構成されています。どの分割においても、10日目のCITEseqサンプルは存在しません。\n\n競技参加者の課題は、**テストセットの入力に対応するラベルを予測する**ことです。提出物の採点を容易にするため、Multiomeのサブセットに対する予測のみを要求します。このサブセットは、Multiomeの行の30%、各行に対して列の15%をサンプリングして作成されています。列のサンプルは行ごとに異なります。CITEseqラベルは全てスコア化されている。\n\n* **evaluation_ids.csv** - 評価対象のテストセットからラベルを特定します。ラベルマトリックスのcell_id / gene_id 識別子から、提出ファイルに必要なrow_idへの結合キーを提供します。\n* **sample_submission.csv** - 正しいフォーマットのサンプル提出ファイル。詳細は、「評価」のページを参照してください。","metadata":{}},{"cell_type":"code","source":"evaluation_ids = pd.read_csv('../input/open-problems-multimodal/evaluation_ids.csv').set_index('row_id')\nevaluation_ids.head()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# EDA\nmetadataの中身を見てみる。細胞の種類はどのようになっている？","metadata":{}},{"cell_type":"code","source":"fig, ax = plt.subplots()\nsns.countplot(data=metadata, x ='cell_type', ax=ax)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-09-14T09:24:39.853996Z","iopub.execute_input":"2022-09-14T09:24:39.854450Z","iopub.status.idle":"2022-09-14T09:24:40.324059Z","shell.execute_reply.started":"2022-09-14T09:24:39.854414Z","shell.execute_reply":"2022-09-14T09:24:40.322918Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<span style=\"color: blue; \">HSC</span>\n――> 骨髄系前駆細胞 -> <span style=\"color: orange; \">EryP</span>  \n　　　|　　　　　　　　　-> <span style=\"color: red; \">MasP</span>  \n　　　|　　　　　　　　　-> <span style=\"color: purple; \">MkP</span>  \n　　　|　　　　　　　　　-> 骨髄芽球 ――> <span style=\"color: green; \">NeuP</span>  \n　　　|　　　　　　　　　　　　　　　――> <span style=\"color: pink; \">MoP</span>  \n 　　　―> リンパ系前駆細胞 -> <span style=\"color: brown; \">BP</span> ","metadata":{}},{"cell_type":"markdown","source":"分化前のHSCが最も多い。細胞分化には時間方向があるので、時間変化を見てみる。まず、各ドナーの細胞数の時間変化。","metadata":{}},{"cell_type":"code","source":"def time_ev_of_num_cells():\n    mlist=['x','o','s','^']\n    clist=['royalblue','orangered','hotpink','limegreen']\n\n    fig, ax = plt.subplots()\n    donors = metadata['donor'].unique()\n    for i in range(len(donors)):\n        df_donor=metadata[metadata['donor']==donors[i]]\n        dnr = 'donor#'+str(donors[i])\n\n        group_by_day = df_donor.groupby('day').sum()\n\n        ax.plot(group_by_day.index, group_by_day['donor'],\n                marker=mlist[i%len(mlist)],ms=13,lw=0.5,\n                color=clist[i%len(clist)],markerfacecolor='None',\n                label=dnr\n               )\n    ax.set_xlabel('day',fontsize=12)\n    ax.set_ylabel('counts',fontsize=12)\n    title = 'time ev. of #cells of each donor'\n    ax.set_title(title)\n    ax.legend(loc=(1.01,0.8))\n    \ntime_ev_of_num_cells()\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-09-15T08:38:53.773293Z","iopub.execute_input":"2022-09-15T08:38:53.773878Z","iopub.status.idle":"2022-09-15T08:38:54.105767Z","shell.execute_reply.started":"2022-09-15T08:38:53.773828Z","shell.execute_reply":"2022-09-15T08:38:54.104135Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"ドナーによって2倍程度の細胞数の差がある。４人中３人で、4日目でピークが来ている。\n\n\n次に各ドナー各細胞種類ごとの時間変化。","metadata":{}},{"cell_type":"code","source":"def time_ev_of_celltype(ctype):\n    mlist=['x','o','s','^']\n    clist=['royalblue','orangered','hotpink','limegreen']\n\n    fig, ax = plt.subplots()\n    donors = metadata['donor'].unique()\n    for i in range(len(donors)):\n        df_donor=metadata[metadata['donor']==donors[i]]\n        dnr = 'donor#'+str(donors[i])\n\n        df_donor_hsc = df_donor[df_donor['cell_type']==ctype]\n        group_by_day = df_donor_hsc.groupby('day').sum()\n\n        ax.plot(group_by_day.index, group_by_day['donor'],\n                marker=mlist[i%len(mlist)],ms=13,lw=0.5,\n                color=clist[i%len(clist)],markerfacecolor='None',\n                label=dnr\n               )\n    ax.set_xlabel('day',fontsize=12)\n    ax.set_ylabel('counts',fontsize=12)\n    title = 'time ev. of # of '+ctype\n    ax.set_title(title)\n    ax.legend(loc=(1.01,0.8))\n    \nctypes = metadata['cell_type'].unique()\n\nfor ctype in ctypes:\n    time_ev_of_celltype(ctype)","metadata":{"execution":{"iopub.status.busy":"2022-09-15T08:41:13.232346Z","iopub.execute_input":"2022-09-15T08:41:13.232895Z","iopub.status.idle":"2022-09-15T08:41:15.830276Z","shell.execute_reply.started":"2022-09-15T08:41:13.232846Z","shell.execute_reply":"2022-09-15T08:41:15.828871Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"* 減少傾向にあるもの：HSC,BP\n* 増加傾向にあるもの：MasP,MoP","metadata":{}},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<br>\n<br>\n<br>\n\n\n# データに触る\n## Multiomeのファイルが大きい\nMultiomeのtrainファイルが大きすぎてメモリ不足になるため、工夫をしなければならない。[上記](###Multiome)で見たように、`df_multi_in`は0のセルが多く疎行列(sparse matrix)である。0以外のセルを扱うことでメモリを節約できる。CSR方式という圧縮方法を試してみる。以下、`while True`にして全部読み込むのにはだいたい30分ほどで処理できる。時間がかかるから`for i in range(3)`で３つだけ読み込んでいる。","metadata":{}},{"cell_type":"code","source":"\ntime_start = time.time()\nfname_multi = '../input/open-problems-multimodal/train_multi_inputs.h5'\nchunksize = 5000\n\n\n# initial state\ncsr = scipy.sparse.csr_matrix([])\nn_start = 0\nn_end   = n_start+chunksize\n#while True:\nfor i in range(3):\n    # reading part of files\n    df_chunk  = pd.read_hdf(fname_multi, start=n_start,stop=n_end).astype(np.float16)\n    # converting from dataframe to csr\n    csr_chunk = scipy.sparse.csr_matrix(df_chunk.to_numpy())\n          \n    # stacking data\n    if (csr).shape[1]==0:\n        csr = csr_chunk\n    else:\n        csr = scipy.sparse.vstack([csr, csr_chunk])\n        \n    # updating state\n    n_start = n_end\n    n_end   = n_start+chunksize\n    \n    # calculation time of process\n    time_now = time.time()\n    time_process = time_now -time_start\n    print(n_start,' time:', round(time_process/60, 2),'min')\n    \n    # check the end of the file\n    if len(df_chunk)<chunksize:\n        print('complete')\n        print(csr.shape)\n        del df_chunk\n        break\n#savename = '/tmp/train_multi_inputs.npz'\n#scipy.sparse.save_npz(savename, csr)","metadata":{"execution":{"iopub.status.busy":"2022-09-03T02:42:20.331913Z","iopub.execute_input":"2022-09-03T02:42:20.332357Z","iopub.status.idle":"2022-09-03T03:06:57.900079Z","shell.execute_reply.started":"2022-09-03T02:42:20.332318Z","shell.execute_reply":"2022-09-03T03:06:57.898168Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"savename = '/tmp/train_multi_inputs.npz'\nscipy.sparse.save_npz(savename, csr)","metadata":{"execution":{"iopub.status.busy":"2022-09-03T03:23:52.576930Z","iopub.execute_input":"2022-09-03T03:23:52.577999Z","iopub.status.idle":"2022-09-03T03:39:30.199731Z","shell.execute_reply.started":"2022-09-03T03:23:52.577955Z","shell.execute_reply":"2022-09-03T03:39:30.197737Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"time_start = time.time()\nsavename = '/tmp/train_multi_inputs.npz'\ntest = scipy.sparse.load_npz(savename)\n# calculation time of process\ntime_now = time.time()\ntime_process = time_now -time_start\nprint(' time:', round(time_process/60, 2),'min')","metadata":{"execution":{"iopub.status.busy":"2022-09-03T03:49:19.229215Z","iopub.execute_input":"2022-09-03T03:49:19.229654Z","iopub.status.idle":"2022-09-03T03:49:53.872550Z","shell.execute_reply.started":"2022-09-03T03:49:19.229619Z","shell.execute_reply":"2022-09-03T03:49:53.871300Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(test)","metadata":{"execution":{"iopub.status.busy":"2022-09-03T04:06:19.756231Z","iopub.execute_input":"2022-09-03T04:06:19.756754Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"上記の方法で疎行列を作成することができるが、kaggleでは他の人が加工したデータを利用することができる。例えば、自分で作成したNotebookの右上の矢印をクリックすると、Dataという項目があり(+Add Data)というボタンをクリックすると他の参加者が公開しているデータをNotebook内で利用できるようになる。本コンペティションでは例えば「Multimodal Single Cell as Sparse Matrices」というFABIEN CROM氏が作成したデータセットなどが非常に使いやすい。試しにこのデータセットのCITEのtrainデータを読み込んでみる。","metadata":{}},{"cell_type":"code","source":"# File names\npath   = \"../input/multimodal-single-cell-as-sparse-matrix/\"\nf_trci = path+'train_cite_inputs_values.sparse.npz'\nf_trct = path+'train_cite_targets_values.sparse.npz'\nf_teci = path+'test_cite_inputs_values.sparse.npz'\n\n# loading train data\ntrci = scipy.sparse.load_npz(f_trci)\ntrct = scipy.sparse.load_npz(f_trct)","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"なんと簡単。","metadata":{}}]}