From bb62ff8c73267cdedabba9909255c52dc2bb02c1 Mon Sep 17 00:00:00 2001 From: Charlotte Rich-Griffin <25774602+crichgriffin@users.noreply.github.com> Date: Thu, 9 Mar 2023 11:20:53 +0000 Subject: [PATCH 1/2] updating with sc_pipeline_muon_dev --- README.md | 2 +- docs/install.md | 2 +- docs/set_up_ssh_keys_for_github.md | 25 ++ .../funcs/__pycache__/__init__.cpython-39.pyc | Bin 0 -> 305 bytes panpipes/funcs/__pycache__/io.cpython-39.pyc | Bin 0 -> 15439 bytes .../funcs/__pycache__/plotting.cpython-39.pyc | Bin 0 -> 8006 bytes .../__pycache__/processing.cpython-39.pyc | Bin 0 -> 11442 bytes .../__pycache__/scmethods.cpython-39.pyc | Bin 0 -> 9646 bytes panpipes/funcs/scmethods.py | 60 ++++- panpipes/pipeline_clustering.py | 152 +++++++---- panpipes/pipeline_clustering/pipeline.yml | 24 +- panpipes/pipeline_preprocess.py | 2 +- python/batch_correct_merge.py | 17 +- python/batch_correct_none.py | 13 +- python/batch_correct_scvi.py | 14 +- python/batch_correct_totalvi.py | 42 +-- python/batch_correct_wnn.py | 8 +- python/collate_mdata.py | 3 +- python/downsample.py | 5 +- python/plot_cluster_umaps.py | 35 +-- python/plot_custom_markers.py | 20 +- python/plot_custom_markers_umap.py | 11 +- python/plot_scanpy_markers.py | 91 ++++--- python/plot_umaps_batch_correct.py | 1 + python/plot_variables_umaps.py | 9 +- python/refmap_scvitools.py | 177 +++++++++++-- python/run_find_markers_multi.py | 245 +++++++++--------- python/run_preprocess_rna.py | 4 +- requirements_minimal.txt | 6 +- 29 files changed, 630 insertions(+), 338 deletions(-) create mode 100644 docs/set_up_ssh_keys_for_github.md create mode 100644 panpipes/funcs/__pycache__/__init__.cpython-39.pyc create mode 100644 panpipes/funcs/__pycache__/io.cpython-39.pyc create mode 100644 panpipes/funcs/__pycache__/plotting.cpython-39.pyc create mode 100644 panpipes/funcs/__pycache__/processing.cpython-39.pyc create mode 100644 panpipes/funcs/__pycache__/scmethods.cpython-39.pyc diff --git a/README.md b/README.md index feb68c80..2ccdfdf6 100644 --- a/README.md +++ b/README.md @@ -23,7 +23,7 @@ Available pipelines: See [installation instrcutions here](https://github.com/DendrouLab/panpipes/blob/main/docs/install.md) - +Oxford BMRC Rescomp users find additional advice in [docs/installation_rescomp](https://github.com/DendrouLab/panpipes/blob/main/docs/installation_rescomp.md) # General principles for running pipelines diff --git a/docs/install.md b/docs/install.md index 0fb1c1ce..1ba617a3 100644 --- a/docs/install.md +++ b/docs/install.md @@ -33,7 +33,7 @@ conda activate pipeline_env we include an environment.yml for a conda environment tested on all the pipelines packaged in this version of Panpipes. ##### Step 2 Download and install this repo - +If you have not already set up SSH keys for github first follow these [instructions](https://github.com/DendrouLab/panpipes/docs/set_up_ssh_keys_for_github.md): ``` git clone https://github.com/DendrouLab/panpipes diff --git a/docs/set_up_ssh_keys_for_github.md b/docs/set_up_ssh_keys_for_github.md new file mode 100644 index 00000000..2161c71c --- /dev/null +++ b/docs/set_up_ssh_keys_for_github.md @@ -0,0 +1,25 @@ + +## Set up SSH key for GitHub +For more advice: https://docs.github.com/en/authentication/connecting-to-github-with-ssh/about-ssh + +After checking for existing keys, if you receive error that ~/.ssh doesn't exist then you don't have one. If there already is one (ie. id_rsa.pub, id_ed25519.pub) then you can either connect it to GitHub or generate new one. +``` +ls -al ~/.ssh #check for existing keys +ssh-keygen -t ed25519 -C "your_email@example.com" #use your GitHub email address +#Enter a file in which to save the key (/c/Users/you/.ssh/id_algorithm):[Press enter] +#Enter passphrase (empty for no passphrase): [Type a passphrase] +eval "$(ssh-agent -s)" #start ssh-agent +ssh-add ~/.ssh/id_ed25519 #add your SSH private key to ssh-agent +clip < ~/.ssh/id_ed25519.pub #copy SSH public key +``` +After copying your SSH public key, go to GitHub --> Settings --> SSH and GPG keys (under Access) --> Add new public SSH key + +To test connection +``` +ssh -T git@github.com +``` +A successful connection should result in +> Hi username! You've successfully authenticated, but GitHub does not provide shell access. + +Activate the environment +``` \ No newline at end of file diff --git a/panpipes/funcs/__pycache__/__init__.cpython-39.pyc b/panpipes/funcs/__pycache__/__init__.cpython-39.pyc new file mode 100644 index 0000000000000000000000000000000000000000..8882564aea2c9024831cb90b1db6fe1942e71330 GIT binary patch literal 305 zcmYk0y=ucS5XU9?BOf#m;-wlALWWW@cj%JMi*RMdQIM@eClk^)>MMEe)K@43dP;}( z4*l=<;}BA_S}hsX%k8EO&CewNOOWj8;D(}QEqAQvHRp*1v8b~XNuP zu&VsxbhYA>u%2=W!E1wPhmLH4x;OFUTO^r-^D)xmwVF*Vq^CnhF$sW>Co6uK6UW{{6HzE2iG;q4g5luk(j_2)eD=}yZ%>Wh*Dj46K{d0L72aQrGJJ_hy VBV4EXQ}5f+ndkHiEGZs%>L1~|O<({3 literal 0 HcmV?d00001 diff --git a/panpipes/funcs/__pycache__/io.cpython-39.pyc b/panpipes/funcs/__pycache__/io.cpython-39.pyc new file mode 100644 index 0000000000000000000000000000000000000000..757b8d63cb336e2331653b42e4313dcc1881eecd GIT binary patch literal 15439 zcmcgzTaX;rS?=5P^v=%C&aPG~tt87W*@`sIW>2vwd-{h>6K77Nk>zXDs zp|`ZQ?&&(^jh5jV{4_n2pO$CwGw0d-%zJr$I-Y~4*_vt>yh3~0o96yjtJt3LX4)mM z)SmTb+hwoZp7Z9~^WJ=h|A4o^{SSHv`MKyF!ZX)8++OmQ0Aq`M>qz^kceH)XJH~yS z)_v{c-tnye{hkZ>sg~6~;hjL;lVV1cUe>*nqR`Ytp<#+yQAYUz@qWMHKe(hJJ%w~y z%r9x;z-y*h5C>n*dG8QUp>7d%hhEdgVX=g|)8e8yChkM;Gvc_o|7Fdq2v?lI^C9sx zpdSGAgMdCI-T~-$iZ%a)e_Wh?%@Aiq1vSg!nPJUC;+?2@m-v8K5v#zsBF>73foWAd zBHoSXS@Ebihv&oMG4UQeAHl31Z@s(ysP||#({tWA@sr}c;tBMAO#GNQ52*Ku9~T$! ze4N+rz21ApPl)%S{Dd$*tW}=PF2@P-~_6OU=lawN}MW zi_))&>ceMmRNK*wbb46YI+spyVe4E?WZB5S5v2~-*28OAZc8-qnp}OW)(pN4GScv- zxq2t4*P?1o)S_Co)eNJw*bMN>(67@F&>{b15L`qnvd-mNv~{Vu%cAzZdWgayje5yHN-WCv zEiv8C1(qoGC-C+);LkN^muA-G9OmO4l(p;ChF^<%(ywmRtL<7On>V~8jb;mr=8JpV z(p%()Ce|gj+kUjw0S#s|YYMMW>vsJ>RKrN7&S2fvg$v(nv75B9pCV!t|?0#>DU0B;IUBGGimuyot7bw{`Kd1JvS7)SKO`Vo|E zw3KxD0H7ua`XDN3-|P5>UpF-npU@K!tk9RVo5n@$Rr94k?OTztXC`38OCZN3keqx3 zWS+y<4!(eT)KR_2fk+!5eS1lhjmR4HNOX}G&TaE!Mn8u>jj=vc=<~&_&lm6NQvlp_ zhP!y%%6gw6yu#YVXts0oUR95y_Vvi#%M8Y`U)94eCpz%_^>D6v%rzf4gIS0J*{ls|W7Z3pHFB8s z6mT9G%$n*(<2ldJeb_F|Ow_To6-=8i4CL^v}xuY2*)QxQdKx z8n}vLDZ^#TCxGh+a8WL!HF0#u0-jGtg}rHUjF=j*EM-_`kiRd>m&Ebmn7)yX*+@A%@GyAjR{^6UAtjGPs9;P(VBP;jI+@mbqVUI_%9`0C=N6{lI_vg6BT(E#Q z&EQ+&H=6x9Jj>$55D!tpiEIw@+-F{#L|sPTVA2P$ucrJ=azH#lrNBY%L9QP{^%p>g z2a<)vPV$M9OeKY6Iw>YINhz64%E??Z-zYz?>Dtw@Fu~MM#TN@MKF|08WNZYM>b5$N zu`CQjce}FzHYr`dR^M`4&A@k?fgAQV+RZS8Xk3YIMBH7HEKQwF|3;N~<$7AA?10(H zMQSx=P-6*FG30%ykjIgDvoUN$E-g}V*f<4=U(fJURu(a+_=OicJ-1#9+)mKiRcjF; zkCg++Do>%)_E|=u3hNm)+y07srQ_B%H>JOcMP7D8-yf4w+_I9+@cV|ZFEwOsv+V~_ z$U{mGll0uE$$CflRfrkctY<9?UB9+dtpyNjHPToZ-YXAbBf+zX?UmK2-Cd2Rst{EW zdF@`D<42{K=2*6-xp1r2^;4rOQmYwa!ydq#k!^MAsUdxczurcghfL|UgD}nEjet2q zJ9d|kovYYop~8p2DUO4DdXlQpTD-1X+UiMFfSkexP-BsGFSK1euDQF4ZohbVa` zCGVnSg_5(BJWR>EDS4EVbCf)W#G4^Pf?lZRn&zQk&{8ttu&h-HiNW;XrXN%lod*G{ zeIWFZUJ`x*iDoV8^AM0WMB_1h%X;ChqPd8glCfy0mLqympEr*{hz?qddKn_rG9E%p ziE{cVw@5fpU({W08>J+6vp<9T0iUfG3; zQ1S>jNx%DsndTZnt?j1{w1{fE*9lSw`%r5&qusg*u)&A&4!(pn{3en}+k;$BG^I4g z<;d8E(wVRd`4KHRW@%Aw_X!>R14<}KQmBx+_M#RTXEjK0v!}^Rk&|dQ^$#mtR$?ez z;O8X!ThY`W#jQDdX#Q>sSfp^&OtQe;1$)lLm-3#k30taf@>TSVpLn>Q} zSJ4&{9cfm{FhHrGWGIA=QNp*No=^=P^)=;lfCO)zBtNj*pG5tH>PU7GeVM!($Y9HL zkjO%!?SWiqi-L7nI0iOxInt?4n1>+N_hE0Kp6sKGS$%G-p2+3}Rj*(nfK6o?@s>O5 zYwBD)i+8HS@l%we5y$WO5qGI#rMc^pOwr4PYsB+V&(cypL-ik^gjVA2C7KXB-%WGPK=?N* zW@>^%O88~8)Q&n;nADXV>^ho@Hv?YV?LibvzX3~dPz-wwDmsu=HKJP232WDUVoZ0h z&6HXjyy0)~f{<8}1}<0x&B>1=Unliq@G+4YQWY&XnM8=<*C8a>Sen>gpM$8+%>4z{f`|etj2mo=mrzbfXR(;a*5Q zB-+AZZV5W$ENodnqDp1)jcG^BVwRy?sSg|*wGdd?LL-!5&h9{Cf?gFyIIS*t1RJQe zKyV{19NW5M0E0RlLei5)vNK+l{?%R+M!2g)G&r8ix*7&f;`Z**K^40I3~Ax{TC3+X zBLUL)p~{2n@Aj?hQLw-}37_)sZI1D#zcWUr0)lrZ4D@2^IpOCu=Cz4XVvP zw1Gnx$#+m7siR~J=a?~o6_DHH69Wk(K)tXk9|r`iN@5*X|7+N0aT}Z_V{Igza?E00 zCFh)Q#1#GtV)`{3WPsDK=j?0Obom9A7sW(h1bd>oy(uxXukRaRPS2CA0g{2GON>Tn zU|!0qDydo}+n@&z<)oO!sQ-g;mIrb)DDlYV24yC;`E^!{HLBn<1Eq1yGMU?S%17Ru0IaW z#7Kw4`FTEawlh7kZv>aX3rH1caV>lY`#NeLLoI2PpW@nLR_mgc)Pv7+ExZB@Tj2A^ zuG55EnNIvDo2}R3YlBtg)+FgsE_5?KB4NELySr|)asje=J6s9fFh5fnn^UyYb>{E(RV0mA&J*%8_A_fpx2-4@4Q8J=KqbXCP+3`lcbqw26zIBB^E2vj^%F#}50jaXe z;{u^)qSxINs*tFq`5ypR{7)?E~zNr6-4rewT>+tn~8Yw5tTbXSGRn+R~@;i9u zz}3G48DX{Htz@YDfGQ(vF9#2dvkzY#w@eWxc;xcrl@Y}yQa1jRsB2Us3grLXb0Wa@ z9qu)gIGGjXB)RQ^D4bWLDJ9x(Yrn%b6vNFXIR(e1GQnp3@_WI zb#)Ugx3{rU@3dFL`Zc;$fs>>YtZuYA8<1$wdZem=`q+t-ADgi$cM^G7*6y5Cucf0m zIt8v$J3|`Qsc`)YJE7u2z0=(#z*AwW%Xn^J*721NPF9#YgDkxDe0KqAQ`Z?`tZse? z4E^JkCk3s-Ck7op0iiIV_YL}~DaU1{+<}v@i3tp|N@2eSB@%K96$ir|$Qv3T?x%1q z0Q5DvIicW^vP=5X0iteVuL*W|2Yl-deVvck8+xkWKvK+Bzq;f5U3K`_z2Tg1)!G}P zcHuNo(ZSQ=6Jfnu9ox#4^YG|kmKTswu<6X;PKRVp7#~(kPUFig=CB$2m_lV%+CZro z3|Vy$U0c1V*&-zYcM#SZjxKyUupMer)^<}Hn+%k&CvXE|VkvNQ*6M_QYVl>8*{*d{ zi*~ItC9e|ZJ|T4m>yuhopfI&}uG3RqA&6f<;@Mm0Fh)O}8f;G{A`^n(g_yl-)PW>nCZ!Ib2$0)9g5k*PIgKc7G$V?V|i`!SjQ zSbRu_Oh+A5`7sz^>X@35#m6)!Y#dW-J16pE$JE9#wGxX3)-S*P!Gz-Y|9diJh&ivx z2mCGEtB5ap+- z_{)^sqU1A3;3pf>l>8N{{VXL+MLd8JwUqC9FOf?L;6yIu05SBFZ#ky?9BQtRi#zbE z;cj~aMpB)Q#Nb2XX&PTLbi%(v@{0b7(SYwm7%w{!Sbt=^f;!Z|_kkLi;&|e6j8mIa zxYkhPyjJMjw+um7V3$MzmpJn%kt@s-Qz+4`ZJvabU;~DP0UZBquG}QyLoGw$bIVAa zZF@Tpk@Xd{IHJVxsN0*`*C3XrlYCPbv$yrV>01y=iBq)`2h>-N`%HbL{m=zKE!ZUi z0OkoTy4&qYV)`4qtUIw+b*MJ6=XXH4FgO6GY-#YGmg>YE*pHwb!t6nZfxs4}(!1`V zFRG|gH>hHAOqVbwJ|y?Niaa78#+U7HdKE7Rga(#FX`F4eLWw5@i(jKqPvZa3 zkA|O>S;8-W^ZQ}H{G(xKE%^sy|LGqMJE`CD2V;M3fW1;+yIpQjzd9v?l7R_)?t7{C z^49=D`V1IGQz$oDJWbT36csl<8f`&8z(1Mj1Az;}nT=w&rb9o92p06b!X7#7BzbNs zMRxZ}VmXqCRTe0@L#?`+@nMtY;F6&U7e@_tB(&d+mJgoVti#Q49qvG=YweNy(9J%|$x5`~X5I{# zZ(@QU7};_#>sdfy6U*^Zu@B~x=!djlHLik9j*dV0ps{jR{wx5YRNCqVBgT^!&Xb-7 z1-^UX+Taw#=dwtX@mUTdf0*Nr^T$E|Ndp;sxM;~u~_|m}}GtF-wtxEGd@dJm~CgNicym3rlgx?+aYj4@`h#b-7mw{(O z<0IlCX8shuVT43EI!Jy-Mr48zteg5?4mpi-`yjo(5!q3G8)jevvrCZtH*f1vVH*xb zxZdHsxJ_hN+@9GhCFr*ghvTrs?njhn_sR)7I};Ka_Ig}_^PFymTwG%_E#bx*A{j!J zoM>=M49IGisSsl2af}uZj~+&l9IXcZG3+@i$8J%p})1P+D5TQ*n)F$?)8z%+3x z06V%jpdF@rgFOvSVBFf!&A}iap`L=9kiuBLi56kv0!4$l4M!C16X&suln*{WhUW&y z-?6)4h~uul!OGaXzQpi0>tQOKX!NcOzaW(+AO9k2QBiFN?1E$j8<=#ju1rtnt8c4C zEi5D19BFczX7aHLg)co0Di32Asomg#9AEPnrJ_vF)9!S(oP7cnF#Q%>R_O+5G9c4(yYU{ zqld6yt8Vm^e~z{ZdLRpq6iGhQzlZ*GrfFyQ*?HO2aBHL~_WsgeBUp=t1R^kTb)cLw zfWY-O-(GGAlWe;4(BTyZ(3B@^)aB1>2ylXJM>#-~CjV(J!wCe98xfay&SbquaOB$j z@(`{x+0NW&mbnj&+2S0I5iyKJxj)GWl zxfVoVZbbLQcHzl)*^L~AwO#T`(JkIbpLzO)W%tR8SC-u?Pd&HnzMqmOuRQtGOGAu& z!+^64cT*3iau{eu_@D|R6bo{baFLtC9bOAIn!Za8l_4H~Wpl-KfA%Eac=8GNqXbRB z%PzUi`1bJROUpxOnhQQBF|-VrQBOekl|zr<)-n`yV=}PyXr#EgTt&PUP@EjOB(QOfxRH*IvsNU{1-~)B{%rLBqloKC!SFl4vS7E7O1R1{NiJ_>GjF$WD zUI9KbISabEqEV6<F%X@4k0Q247sKriR67phAkD6W zsIkWQ0rpz{yHx8+ zv%}jQ7e=-PM?%E2&yWK8^bM>ZNA4-b{1LFfE3BxmsW6EU74q|xkj@^LpA4erMn~+r zPr<*2r@ivYw$r<7N=mGIV{qaT&^QhptulWS$D$_!w(tbZmyA zpVoxPiM+K{SPj2NFMop)QuScoj7%Al7l=oQ{C21=HpmX}4h~IxB^jW%TQjK*$*>HTGL2+^NGgiBJ*wm%s@@ znQyl_c4ERDB0G=_LY8b_#tPjzt(9?OwA6v~=YkvWM; z{C+RLui-ZeW_F!nMg{IMO}u>HKrl)1Ch4s+pOHQ9X1eT*n)qJwj4r zZ&t~yk03pzUQVx%fOEgK8T3(RY7ZysJ#$cya* z!!g--H~snt{^yuGk}kfj|0b+;Gx7UBOoM8N^@e-&@*@T@%x**t>!(D6sxF#AK*X&*6f&-^gUSEp~>RN�uX$$ z5zmauh!9AP>l^WG7OE3cED$Q?3q%j5d=ZZb+l9=$>q;9@h0w)BQ-%~5+D;42kW_RZztNDtgl?&c zlKd;ezC=a+N@{iRTZq)6OZT+EluJ>?@a{{$-NEh99}=7cEg!+J(EUF~eQ#01E{rcw z?zaiV60*~+ILc@^j2bV;Lqu$(l%CbVIJ{hger4i0gQb!;sQZaw!WKqZg&xswium^< z;kW3>J_I2ASc2oFvC z3QC)lSD+oj;(+$_NTyNZ(lRn`8uB67j&MLAXp4UF1Ncp?VIe#P_FtyKl@q;*Yan@< z!bix}z+oiDx+EhOhSgu-L0dq9Lbafd>lM;i(gP#Cg?@rU!449JwM9%hbpro-m2>(| zs$aWnapXS(rEXw77Lm~q+e$2q*bl7<{e;BQ8D?VcP|Zem7HylF^PbN0UdFu1Hbx{? z1U`W{Ea=dr^<6E3OHs&AQJ_bHkhQDG;Y1kgskII@3l_s?fQ{!V|A}(nK$0#%e76v; zz>^C$s{uem#H;B-{tE$7T*Te;`Zv_g)4U4EP_UAUn2=OIF&ElL^5ioo{U$QxSh}f4 z6va57z(NTAl(`It$l`|^bXlx}`>pqlxWw1O@T2#2gwC!Y^ugYCb{jglEL)-ND@r38XIj^il*9qko!jeka<=pN2~ zj_a%#gQt8RK=LC5sk1W{6tq^+*DKZ)W`U~kV_-SZzdLP-7i;oMl$@vJ5E2?Y!Zbs6 zbf{Q~PjNLC2-=Gz5|>09o*m++0Nq`65k)6nxntQA%RX76FwdeJZxuI6>J0TFWls4` zrYNQ~>fjPd5~-yHcmnIBrQt<0gpvGYUOsaQvquJj_q{gyzf8?^2f|_c92_oRK{2&9 zfR6)pIhvG%C^-;M1^e(F0N*!!jyjVbr4%M|nkcP_^a})R)|=g3XqAeIe~ns+*YTm} z&6X(@Xh21AsD{W6$3HcET9u}eH1%!A(NJwaPfaNMJ1LBDXz$+8J|5W#Q z#vwh^r_SZCQ|Fxe|Nr~`zrAv~sNwU|PoA`2e@WB+nHr}*1{$~UP5uc5*EkEc22)qP zVW_Lww;DFf$DD?XK0Pe-i;be{Gs2mEsZpAYmm6huPo+`8)eL9*)kbwPHrJTrC2l zD{1l$YHrt0o9(_gK;0N7e#3aX+kNc3c&+p1o2}QS*0$-M$;ZXDgzrz_n|u#Nre|yo zC1*EOxu%V^BW=u%w2__JV?EQ3uw@#@h>_V`Khn=go7^~Jr+a!<&(3hnU7efUIy{l? zIh=8uJBRv6&+u}Fn)q7k_6l@}bbN=#ujnIZWTeI33|6K0O1wlX!?nDISDTrY&P6M7~^Q-aNn#&UL@VNm-^jU?7fSla$>=K z;qyCnUFxkpS@4p!AMqgS*7al3D=~vYvK`5`9~si=1l=Uq_hq5uwSCMP3t8}aH?Il? z*<&j9ofo+wulY6Nr*{)=J(%2gM(iwQ6pQgTKiL|_Mu;s}@7)1U~ zz|$>hhkh69k-E1ljkwdP&q{k31@~~HoyK91q|z3CKfaHNgg*$eDrtL3x;OBpfsM*a zLVA(!IPeoWlf)uz-uJ>G8pa?@rJD>}bVO(-5eYy#^#yfF2U{HW(QA6UK_Uw@SsJ9F zUw37xE8^jx+1`qSw%;g1a#F{Z=)`$P8s2UqZ91l*ke0f!QI69sUra>tOq}Q=O}vg{ zdVF`IJLn{@ZS44AxY71R8h3>^*xK8GK8s{yKk2>k<8N#vqP;QjqCqh5ljdO93WKNZ zjm|J?CmZzq6z}lLVDDr<%{cOB@sy;7LSrSp%8L3Tn`h_site&S_8haB`MEqlgW&yT569Sx}&QfF<#s@#cyOz&*d|` zbjU_dWM@tn*NHKG%*$OyGmi{j8M$cztHINIFsE^eR&-nY;31w;7!@(1#&{?;qjEYsuBLNY@qqE!6%DoW0b5ayA)U|4FxHvO$y^w#krmcx z{257JJ%kU=TwZ)r8@czsfA0^u1{0mzY}i9MuI9DL_0?lKcx8Bkyp$+YK`)^Djo+rv zpWlSF+FQv$Y5R`4ig}Dnb_X_XDVr7+-9^Uum}PZcT%+czDC!H+<{}kOYHkJydUG<9CAZferr&}?XR?ADibT8v&k(_i67d@D z5Lf9oBjE5My+%2|eQJiev*VptB@4tAk|xW_iyx!=megaK6C-b?J29Mrco{><=TK;@ z%3PR6mA%9&Fpi3DGZ%(pvhy&J68hKKB_@6Wl5KL2(+|Bd`hE#tH4!{A421k}%6h<$ zj6uoa&jXej-)3J$KR595&HdZAg||~X<#cL0LAq5N#6gtS$aRRM))8^P=GDmOrarIT zh6}wTyuMH46Xj_^NpYR^8ZBnS|hw7wd$|sO6WLl5>0AQ1R12+y0=)2a_iANmLmPy~O&Ep1S z;TQb~Is{dAMVxwR5Jz%hD;B|i9Hn09g+Vv!Lj)Ra?}lE>53#jb^1!=k(B4)(iJT?e zKiy$KqRzwxJmOQE%Qd>W+azy(8xn~kR&aPa)8YEWucEI_kTv~Koen^4F+v3n51Z_< z`>@F9Rca9@nDo3t4sk^rv1nC|7$+mQ;ONwdI--d?92-Y;Jm3dw+O{rwsWrBN>1;V? z8$@-SO=*b8+c!~Mm*xN`wQo^*U77a&`pZczNNQdf);7J(OSQHaMR8gq0H`}a%vK_P z6swkIkdVD6!@jhVEpOn9cc|O)g^1}DSbf6fh67t5s_^W|t--$h_BIVQNg117wAb{K z#M@H}Pre6`f>W$QQL4JDpNI8`0KIKekm-j`yMV7c%?{KjrSmiobb<-6ER*boRObO8 z1a8#q7`F*TEH$cF4Z%lh;V$by$5o>%_>jR(l;BxCoFLd08`Oqajb9Z z3OWc_3cX@xgHgG>45~WAXAcQ{oZbvZO0)u8%S2u=Hh*aF1-^J_L8YBur8f&8`FP}J zZm*iT<2k7O5?>~`LeD&UK#rgGj2DlzkC-T7Yzbp%jX`+;(d}fWFZAkCY@OHr4HkJK)n}~aU@IAfK?<(x?}AX9bdgIW z1JVxKoPe>$NuYeairADuegU=D4)w8>Ax>SLd%P|3dF=$ zXp1E<#<(r*JbV(j&?g;)Y6Pd^cc_JM1`!qFx15|B{U8!=(?B_C$HH$WULRJBu*csW zh%aCOft0?xE53!Q_+1iYzbuVy#H7Sm8uMU}O`f0UG-z`s?PuR{v5DK7?V-5u|0eED zh&pTRB1qgNyf47f@1w1J!5K^~C571%`RLh_}6)zL1X?I0|2xO)!hruvEsiwap0 zEU^Ien&E6zN@p^gf{~Bd2d^Of048z5t1>UHXnZC^X5i3*lv=^4g0+@2XHlD!&I;53 z>NA^`vso3=^bqz`M{@{q<}#=SFRNU`{AeL7^=7kqUda}4cNKF}c4027u4+I%!aFys z;whMg);n5+eD$<^R?QYC@(Uz??|YD#>&*i-7m4;|rSSq^%4gs%LMB4VTmmRKkyylO zPGpW0Z^m8A*(@;vk{z!cX+YfPMoZaJZ*{!JSJepB*ZCUNKVtX3dGBl4%=jE3AA+%x zQtkEqFBfb0ya=+){aVXQz@X!(ww0!XSS8>;6a56@ryLYITX^)mxXnxzCS`AT$Mu<`r6M$$DiOsiJ*e=?V;aEkBOy8yM+`& z+^?U{0mnczlH0r4_JdHorDOnIJnl;~@qI1~bSYcULk}B?jiKxzeIeHl3^m9cU=zL#P@iTDF- zp;1scM_W*DBlllFOPl){Dt-+Ga%4fbnWwReRLyBd9>Y?2B3ab1qb5z2q&h+3ERl8R>+(Kg`(R{?UX zh;Qdv4G?ybKu#y1BL`Vx^GF1Vzr+1!AdB>qV2fUG5~BD}`%s6UBA)I()O)Z?cm-qs z^$%Hx5fgvl0+J8}Qa{Bu1Xz?-;Wnb2JMcWJFEigmIB06Av4>EA{s43irIz0Zgx%Gy zsnT_*1OP}TKlIE6cqaHUmzl*?D!O)J&`*@vFWTRZdA*xp`s}EDm*Yl z-WLHe^1kLnwb-0$*V()IS`D=_Wz~}1hh^I8iIYg37AN{B{){T0y{Z$It*EqeCi`XV zvqsTjhbVX565q!h(h5?)pX3Hz{3#-oKS?p6M3f1ehtxmPlZcRseNMWTK=D4Ael#@v1Rb4bts#Qd?S_O7`>lTI0v2NEE2MXE?7`T?*9u= zC5==mWl*I{P%}M6m`E^v0A>W}ajAu|TQrs-W2pM5ZcgiWHPOSQbrbePd15n9BYT4Z zWtuG*V1j*v>LB@Lfa-we`vcXFp$JFRWzK(!HANPOH%fO5>X*XVwA2(q*4ub zrIM+xLUwCp@0alF1+oky`#mjMg{_%;S#E2Yjijmlz82|MwbV-O6e-or_$n+`afPoy zQyrB)MN0Dkcul`fjPd(~T8&HVJfeI?0UX7JWtVW|f1MS0c?CMnD{Gn>EA`5~3b80! z8A^&{)hLpnvwRj@#=#z{hc571)yK10LA|Drz))<3-*ca~e=@AQ_$C)oD8pJGF&r9} zo-!B~k}naopgW0Ta}CwSzd>63I~ChhoDF1^QT;9TE2DaW{OL3emYg#8r-ntI`Ham| zWOlMNiZRvFO4Fs~cAhVWDSP|Ob^QHw9Z#)Sdw6>h8c_D4Icz z0K8xeVVapBZLC7GdR3Xi(`2)X0A%IuG!;Q>nEJ}|CxP0>HQ(%<52t(`c(`%4Ij%VkV@~YF_WKEQ#Hz`8a{8h1kXHx8+?ssQQ)ix zyLk*!z>h?e5cjOg$nqH_`5WAtTtcA{mf}a1tGW&^K!5AN*Hv7RTZYr4KX`~A+D_FK zPat16HuG#f0b8BnhfVy8Z%;D~b13`74dxJxa$f&$?>&Xm4Xsyg_Ks~4AxtFD8;`QqggE4syNZpoebKiQ89ZvX%Q literal 0 HcmV?d00001 diff --git a/panpipes/funcs/__pycache__/processing.cpython-39.pyc b/panpipes/funcs/__pycache__/processing.cpython-39.pyc new file mode 100644 index 0000000000000000000000000000000000000000..cda0969931630816365f844b8edf02fe8b1aa752 GIT binary patch literal 11442 zcmb_i+m9R9d7l{$$su>OUai-QPzqb%c|`(Tc?tfnz7By@(g#$ zC5P0Rp;pUsbFq|K0je}ivo=b6m5Z`&qdLP`~?L%uLTPDr7s2gU?6ROzweOT zrB(@ypxm7^=YIVzzwbK}95_(4@O$&;Kd*lK^Op7ROzi!+NW6?E{;6$Q%2IaIYT0Gm zX548yWk=qbaz@^6*_C&;oRxR3oRfFH?8&=OF3NkVJdJmzS!f+7ACSIW#AnKnNj%$} zZ5=EhY|WMDB%NzM-kLAZOFZ9vqIIZzNa9{|s&%-0SmK4|k=D`jQHd9u$68O8pOpAi z^LXp2@>7UUs{?nf@(DGg9>e=-HLDKd?W>o=!{IY(?juJ%uIBG%%O{b3Hhfk+@li${ zQit!ljykqxt0$FJbJ6PfN49!Ooj|M8>S^VpzccC?brSEh>RELP z?~-~>oyPll^<`C3&!dmes8`iFwSd?;wWyX*wxF)63+e^rEUMSkMfD%`u-Le$zKC&`)k~P$i)+>zPNJ4kSJao#?z76dZk4Wm%rh+6Y4Pni z)YlVTkJeU@%vIY_HAoQ4wnDuYmYme9wp$$?#&Mdht_Ssqwc4w3Mb!{1#NE|6Xmy$) zv%lKC8YBT?`B$Q77|rb#>ri#8VHI>T{xR<1c@(9btOOV=&eO6;wZ zw#C19u#W52dC&z`;PzfA`uN@K)~g$S62^(2wEZ{?{d&z0{8*OcSA)o34gI96qfq(v zC}{_Nt=SF|X?x?^3Kp1^jI{XmpxF&y)4Hv>7Sq84OV46fpS-=a)~Ur8mo~#@bEz8W zq`jtt&id9;7uy#v_2R~*&tF=Kb#dTYu2`%SPgsN^sOCx+vz*X zK8h{eaI}B3eA7kR8Kz(E+k=d&cl|Dz(EXjFMP72p88Dee6!MzLd0=s_S7GhyF`f z{BNFN-kEd$nWW9ol$?3bmsa&C?$kjy{_2*$o+O?4^3qbA1l5i9tx#jpQsXn_=T?ALU=<=5)XFhY*sjZ~=p zWIgn2?Iz9%1n*bty4r2w_#&{tW$EXQmABr$;n&#YvJWA1e(_SJ+U&;YaWT1_EUxwZ zw)U6D!`vc|Ua2exQPA9q!9B@VC%m%EVg0yUtJQB`+0*a1S$x?)iQZ32*>Ajk1A`(5 zgPELHyIO-#LC+v&)eht6OyY00^#)}L*B~8XPBJ`6d$sM?V?PRSi`0>2?oMIRKP%L3 z)-zm3Nyh%fem|XLXa}d8!*iWSeFVW8nD@O|GBu7IY`Nq7*tzLq`|Z@OKwLRG)Q3>- z0iuiNQ~UPAL|0sH*5iZ(euYxvzwlVS^O!16tTNt{3J(n?Y+;1y`=2k|s0s7C-HTOhaMI z(v*G*ZPN_2jxbH?a=x1TmmZMy)Nx9|YBNl8pb5m&ybe3fpc;xO)d!?Gq(CVTxlgjw zr|7f?!F)Kbu(soQF+xXK+L{Rp$h5}}X8J*xrTltAt)dk5;WZ!4#Zy+wXz z;OWc!c=DzjnfOG!?-ckB=BB%gSV!XCLGvg3}k?XcZd>%Ak8BgMwlV4CTE2ATtM z0-ECM_pHHz9VlUe5Bc7I3prnh_74*uaDcoaBG^A`6n;#`twLSKC1`~OIRi>@>^-DK zST!0r2{tGhCORP405|&WhzuF|*oy#KrXI#VPh0djGvjj{0>X7~4KOT=$ek;Y?~5|c zx5iRRk_+RXBxTDl`v$y74siv4CJVD`jOP2vt%*DivkGk#szv{e#BbHt))V3l9Rik6 zZ9!Y=a6yLwqL3q|Cx5dJRrxK`@TF^DyzpDB*$OrQRZwlbM^w0!LyIdPlNfqufoGz2 z;_Gg-Td9KG!)D<>9PPx`1oyHj=|{MiqC-T?%(e|KhVE}zM!*aiIV}*~ku_p{0+qqq zFW|My&H5Q6LuLshpB>MZ6RE6jRd4{0k|?E_)T_s4)v31`=m8&!ikQ`Ktn^j*Bokca3gegBR}xarqCyP%a8}bfkoj!DuZnBS0i2g zD=;CMR@G~@5NfX)3eA5SP9$v#w|w*Vdwvips(WxVIoM=xy+=A0#)Tn~z*-oRJs2|b zZ_7Y1H~`BLgtx_TSoGg%cm2%(9FMl9M`;B)JErGxxAh3{K#`MTAyET6bQGNPSG$S7 zg|1ZFCy%g4v$Aoos5b6sG(U92xp8S2P#XOIiHq9?<>(hd0KYOz>Pu)XtY`>VNSOLY z81dVT2qT`xsq#?wQ1YYLa0OcQ3v6_m!50~j&-9lVyv*Pg247}C{9yRzDqbyJHl1>>}Bgi>qEN+O{ko^VH`J4sGhge%WMgzGU1*U>}< zcHdq@mq+&sS}6MrXfJi^Fa}GGei}{7S;4+QcG*C1KpwVurP^*jBz%8}y5`JJsrF`M z6mCVDJnF=WD?|jQ(;8e0)CUZnOd8@2y5u=r#}n_>5irN?25`h49f8@FI0GUYLzhAV z5{&8c4_V0lC>8Y^c%>PvB|W0U7EHHFH-h1;D#!}g+8U}Il6oB1l$V+?Q;?m0o-vM? zy47}PtCZKj&V&<1`m2%)t6!x#SwGYVcCMV?2*XZgwQ4r*E9`(XOJ^V6ukRz*pfX~y z%%RK)yo>g1W}b%GZ?R1VV&ew1AQeiLM?Ly6&WBtQbIyHN*`(FJHL!rd;JWBTcGLbQ zjsm@9L2qu~YB(l$*-C)1fYDsgGq6KhnBzs-QerppLQ zSd1`kU|lLwa`+_OEcE=lmu&@iTaumr1~gxLzW z#JfmOZ~+a3EFRekl23xgCk~cJg@$dQ`FY*a&!PxoM97Er(<~@*E!3Pf44lGhyWOk^ zrhjN=h0XPPb-hxL)IXtWOmeksT|qbjlSgDM^O!*?=Yfr%7=K8<93Q0oK3=H@Jpaqa zMe)C&-e}wa@rm_w+zcKW_gjc4am(~E`W;7qxes3!%<1FG*|u-Ki&O@w=atJ;;ta4$ zgz$;0vVH3icI*Q?Iz~TOl11ycNDL^;jmpF^lV3^Q8@e0bNV9? z&``(5ov>Q3ZTSle&&7Ta`)onfEKPcE?6&YfT@(K-IthFplYsA`pRH8Nr9}!dW8A|) zOJI35+iPpU1`u&kt(_KbV0z+SOI@}vIpt~4VmE;_B7}fA>?Vd^$C4zSfjg+WQK_|c zg>8O@ofRY50xm0xlMIeBr)bY0MaPQXK#hGTGiFjkz~1c0!>H_S`x{nt+ISZeyD>Zo zoD9w>^04P%&JXN&tjIYJF$k7)SM+)Gk-ESsVee2R&_Z_9xCUx%I6lCXv@n;V z+bCuOLA8*^F5Q*7!QzL=e2yg*)~#U4(&xXSno*uM-Gat3m~&=q-9`yKu!5p0b#-k* zk&?Q{9Qjiugq}c+@mvB$DaIhylw=xk4aKLC2M+xJAQa*$+sNIr^_3*AGJVI;BgL2G z-asy#io;yi_1HpXmqD*c6<~P>b!KA|5un&jOuN5TK zb>l!>4U^3IEHIu^?*tf& zt{}7pj2sUGzhtx9isIhkiE>pN6_kssfv&byh?9dGBdjGv1n^tJ=|~@Ye%QfEX<92B zpQdM6uOoo*qte1v?8-InifC8eWr6pVaLaI57J@Rot*OydSLt>qbqSY2nMQT0S~&-^ zh|uDjtdonvO?dPlCdYZX1gOwnHSmJi@;}fZ_7PZFvLMVNQ3^$S8bmoK%;Vc4 z`B4h3|CRL`Q^=tvh9~vDOJb$1HJBdt(KtZg>3sk<{WGbb1%f287JwbmzLHI529Gf| znB5Tp9?c{N2Xo2eDo+H->m`u0x6*qGUVo5>p^|+Rrf*^cv_??MKxl#zD(xusG@S{E z=@=z9O1AFd-Ft~$E#3+=y9?k6sK*23d+ajH!s|Kel#8zXotO{4V5c%tH zdF)yk=@hS-zF|y$=n9vs&~#(z?=tuvg0ja<$UL1E0SVWJy*?x$rvjx|r_xFkyirzs z3QaB11+>qn8FS77tH5smpIr91$V)ic9sRqgw+}0fbv7{zjuJlC9Sc_}I6fS-9H7B0 z8GQ#e^zR^m$Q`X2@PYhPK0x-7UD!qIyGm?!=$zoRrP6Uudy+%+Lmx6!Ugy?K$xV*wqA z0}(<{1gd~gq4g;<8g4G&&7g=ym9KNc(s_Yh$wfn&mvP`2 z{zFJOc#Ke+cR#jOw;;>t%U%6}Tq4{PSp@qSc7}MG452_}OCMR1jSUML;^90S%Oh!?n`U+7{CO2RF6R(&!|~2qokO79joTN02U}+ z5B>rICHq{*GeMvreb}<}8%)}>8Do%$gYYFZ8k*y&-FiqNipS;RJ}MRG*CVNPALDRd z!Z;@p*>x3=btrmnSsES+GSQ^HJ)WBI5a59m50&+pc`lEpwQo!^tLy~xSpS&82?oMTW4S4ty*wI2 z_+9@2s`ZMO$C!6wt*9#8BWopGq8-lL5zvXV!mk$S$H(U7f|!<%KE**N6^uTEyf;UK z_MRT{^!^hbUQ0cR>iZ@^xh<36eV)QbC%C*FU^*aqAwEv&fwS_7?Zb!sN2sJfU_f#4 z3n%rJ(WD;OsTL} zD*h=ezQBsJxTC3zk8A`p!DQ0cN31s?K+R^K$9m&)D5Q3uf!Dq008}Ex-5JE`NT5z< zy_^r=dCED2>va9ss6V0QU1ac~4=8>fLEqMgr~_+u9Zs3Qw0{WU8Q~&rGN+GkOpN9? zeaiaaxp;scQV-6nN~;T$mfOJX#ajGGTkhaqeIR@Cn_QBI9%#PK;o?vK4dN5zLsY=& zzKTaIT~3vE|8lEOG7@lt$qi5g)U!)Cq$6zVdA!mrZjDs@5N#8jscz>T&VplD*hu}i zcu#B?_0Q-L1tyw7JgY|E@Ea3wRA#>)l$}0~M=di*=2h063@#NOZC0(kkD1AI+5 zy6-F~ujMRWN5nvs2mgQw^gn%^^CPef$OyZMh(q!( zi_;H&&@~(WOkjg~alLJ4qvuT;A)G3(!OEs4QgA$=yqIJMIee935Dd`AAio2o^4^Je zo$dU#N87Fsy$_qsqu0Kl>*ukb9e%3Bt=@%-+K3xi4?_DwW zb;vKl;(hx~3#t|dC=79{ydK=5$rOh$(ARN^%PZz)uq7Y4z;_l0b=-UnYfHB)z0b(q zZH)Hg>t|YA5PG4H>&!LW_rhynI>D(BWTaV18ei^TAmfp=Cu(G`7FEz`G*5^Y3`3Rd zPwWMBO8&=pN^quihMWbT`{?dq_voKNweOAD>O`wsf$d1U{!vT?=fg7Un+dzBUWHBA zt9Oik5p}aBqvErRR;&qMHD!mP)SURQVh}UcK)Seg}gr@lc_K8`w;3mcGn0$*3m+OxQ|b@7?i5fJ#it^XP2B^#gotq7RF zCUuZ=4e#&ZEdj*9?@2B5x|G5+fO?Zp-(pSS48wRUHkPfdpP0&gm6DCRtuCe9Zw>VQ$rAbI4!id`F3MB`xril?2YsSX0=Tk6DT_Zu&c0i#QtD)pB8&(dn3r z5i0XEkFNpfi4xyXnlqO>dWuB_er{-boo26d4309mz@Wt78Uq3VxrQ@$2=6jUx4I4( z(48uNM)Bo{4XtUhiaHc!D=JL1h{#dnOQ0B$?cQ@{V+u^YWgK(28D7DuSuss1`j966 zV94|PIEMIRD&tLgPk6`Xo}4>ceB#(?&&A(C??CZb@wn#}XNxDjyqEQ|v)=y!rDmr3 literal 0 HcmV?d00001 diff --git a/panpipes/funcs/__pycache__/scmethods.cpython-39.pyc b/panpipes/funcs/__pycache__/scmethods.cpython-39.pyc new file mode 100644 index 0000000000000000000000000000000000000000..e952e5a5c46b14dfc2b7f4335c52dd4914734c9a GIT binary patch literal 9646 zcmcIqS&SUVdG6|KdX7D~T$0PfQd5>mjzum>>t-mCvScc9q>T{~4$+p<>YeJD?b)90 zVRaApaL2oU5!AQB(}g1iKhmG7^b znLQ*K$U}NCHFf;;SJhwt-4jer)igX`{ph3CGf!#SzfxoGW1(>!-|QCCG_G+LYaWxo zx@X{5kIkOtS&YVv*zQ$4N7l`_+N*grS!Z$GtIIhJuYq4Hp6X3|)4dsQW;8zQ&C2mZ z-W-1I_)zb#cUX>9;`!bY@5l)5sCN|OPJFDl;4Pr9%4_lQ-U;u7tk>g{y;I&PS#QLj z>pkK<(tFf<6zwTKeOL2bKEr46dzv2#ALDadtgi9He7n)WsyeF{pU20wBPx7ZQ?>uL3YRwlO63v&ItaPF*>jxqWQMc3O%;#;? zW|CC&f;id?{XETsxT%-UulMsPO@g?rydLGDKqqFrF?c1&1Jo)nCCRv&1$V-B)Om;- z)G$;P5$}ZkIL*@O7Y50=o}~E&40S^8avVs`&1vqPP-Hmv7_L1?T0{nSEfuSo%OerC@_55V4`ElN zbbCQAqIGw9kh?&U#p&9j#3+dOn^vhO{nEg4r4_FCdzVX-wo<39d>W7)%A1uewZu=tb zHoKU$dSSkj@@%ocQO;mRpO*90UXKwmbrhN{7EpeI;uUmCil0M;D2}g=uTg;3oGqX% zw64k7aqT3;m2N=N|F8ak=`X_D0}wxR1DA(M7LI6s-px|C7j;&0cNu+QFBKc^${^{8 zjcam4U#2@z*D~mb!fgjxzV~j#42~&gQItkD=#>WORa)6f0N{BVqPtNAk*=Jc&{J&# zeTtkFaSHX$4}U<L!dH^{XJGSO14%(#$g{qxZ?U1Z zr^k>zHX!uE;wCrW)!RBCtyQyYaEsgbS#B2A7CWg8D-z27OekAIIk~lM=ap@zs1#0H z-(lO;9bDtET2#9=ULlOb8tU~?y{n=RJ> z#O#57i52~OuyhX(=Zd+)z^)JT8n55g^7-O$FT1ALOA6IcWiP&>`6b1|GHBp$L;!=(a`Lb4hORNF_^(F?VyEb}?P< zUJtS?*vKCLjir~aPgL;J6TAYXgSv*aLFlDMj>9|QXqUTds1rT2=#mX`vpg7`@!824 zmxI=7N5J%OcR9$Q@1eSJ7{A{l@}q!tQ0JGFULic#<7i>18L=7)bsNtujxRwjM1Z%rz#f3s!%8gvB?Sf;w>S*AlF@^mfhgIJ`tG&+vnm|LE`b~VZM^JR2 z)!sh(mN8^QeP|RcT_|*Xjo;I6n{Vzh4o$Pv7cZ48DA{t!S_iF0?P?rlIf=}RJz)JW ze6`KfO89-H3l^`&X)B1c7f~6*ZaxZ@qEn2SO0rx#%1B6w^P08NB&wHGz&`uOQML9(&~#S$;Dy9FdOvy@rED7 zY0~lKRgz6DD|srfu{46^OnNP){?g^r*o^vR^-duC1WRF0r6pO?s|BLt12@j*O^@ny zAZkxh8?IQpMLk{>&Xq5bdxyrjfQ9^XO* z9Hw>QAm7mh9KfNGe7Rr+>*^Gfbj=;O%|lZ{Sl}kcb$pEb=<&%V#q8tSHe6yEx1eM&#|NUm`p_IdH5w1Eg9bt03MR=NfiCpXi#V&sbt(YhI-Blotwwu z0*iN4t*fh`qH7Exk=s^bAnv2|OXXbyVFxtAyX4wKj%gl^ybBFLj(pc1TFLeOSpKg# zh8FIuQULSN&Yfg7&3PCea+c_d`BBftCbH z6xbhO-HeQ(26f(m3O~gb^jVbiy7RGan9TXuHpK?U_leT@@=3h1w@}=JIv+6c`aKB7 zJqXS=+kr5=r@hDC)7u6K&RwTa5fLIyWYY8=$d1U%%&h`ZV*xSR!e?kFjl6PCyQd>s z)Q5Iq<4P*QAHWZj9^b2B?yf+7;G^vNf-Y}|0a7y)2QneI8NPyxcmR17<=Va0cFAAw zyUQCgo6-+k?TEZa5r~^cj~tcWgnPvurMgDDltc2qF$5he+wSC9N*5`lm~0Hru38t-2D~q@*yxN+0$zOMYy2KF_ z4-6$3&iI59TJbsMQN&pke@higGi@%C_>xRa;s$?Az)$mw$F|;TtprgbHh}_`j`VVA zgGFu+!m>&^%N*IoOblt6j0TM?%*!d_5?`h?m*Y@wA4*-CNji}rucx%-zB6xnBB;I; zdcx&ZC!Oj*-=t|WJ9Of+&+$_LDXl%@uxV)Q24uIv#1{JZ$t>jz$^1!~{byA6DQw7U zn~`MRb)dSDaoJ)re_`e@Y222I{4G^OlCA=^$+{MlsvJ9@tnK9S{g9sq2ic~~6Dj0M zUg=o2nUF-izcQa zFqHSW8_FDD#D4;iV#G^yfl_tv0{d+dNdN6qS#iIe5q+S#2}$@ZLg}Lz8Je;h!hJ&; z9c=kBjU6yNg8B~FnHp_I<>&?3y0RjQ0C3{ZFuSrjre;tutZhbGO7 z=g^5L)h~#o>hSMUO)~ANUv}ny1X%Vs3K$;55M+AD@XQ+%vtu}{CQZ^b6W<2h2Y>%} z^!)p^hjhf7`@)MldPwgO1+sH67BK7#{DZ9Ewg&#u;J^_--6Yzt!kT2LZ6Isw&*Cr{tE8hj(*Y~Vc6qRmeBvFuyZgm?$%np-|P@RIB z`wY-^xSJ?}6S+s!8nWF$TkAIX)Od6oZ;+&e+Tb(jlT<-&_3;9!2##-acIgI!8!HYr zLIGEEkjD|MmM?ynXpDC%73sOAiR8I`QUoeXUGgT;oD$vS{?YmzCtgJsqiHP3l!^Eb z1{HJ6Piiv#E!o7m8V9{)9=!M=4NA5BgiOLj?F~Q5QgU_u@y7OwikhXCIE}fbGlnWF z;W}Jl9(oI71fz}#BJmDpcqfv{3dm4Af5{&?R#6-_r%SVy_BX^miTd^+NOE+G-$7Fe zloTCHrWQ&wLY74e1&l#1${K<$k!KXx$s~qy?nBIN)=RT5(wyFy$+Kl;V~pHK8qAxS z>+ck>Dkd0Ul8+1 z2fAr1-?$zK5p4L#>ISj$u_^Z>h+=|bwI)-lNgiwWY|brQppV3)4UFoAi7vJQUmC4g z2$yE+rNOghu_NosbthX0WqK1oWetEna+aMa+4nZZ!r(-wY@?ywomBdFEIYW z&wtW?q4x9txc0#}=ih(zUpDh=Es~k>Lzas;Ob*GP!Htz>IJ<4 zyPB@Su3CdW42u+{=Ap7iX&u-Hyvh=!V8owLA^p8Sqnh;Ku2Jn-Dt7Y&cAm;$fp|qY zRLGaizmFPhRgA2)lczr1!Ttt>lnQsB<;cfI`tDbe6 zf@`nFkv{^ye#k!|3KLr;0YJgCBnnr6OoP9s%t!+!gKB4BHt-9vkRk&A4j!Ivvj%z^ zXf^zK`NY z^FnFk&kQ`sBngp5A(|*d)5uy;e`8U-tC6&Zv6h)$3Ko?kOD13ZCKaULN@pz)^g1+y zH-!IopjO)OO7JU4bOp&X%>G7(KVS^tD*Y*S(#;~{lgLDbysc7Cv+@xB(?Ne>LEaaC zqoA6o5~M>_dgCWZBb0iSiUw6`V6_&#sl><>I^_Ei`L06#FQByAaf*M5yh=lq{FK=` z6+3TH^V?LA<(5iQY9vYRJ@>lDnqCFF4C3%bT&b4bKtZ=h*?*b;mrentV5!xR{j;W1 hb?OUm)ZE%sZN7FKpOa2iWB7ci arg_mindiffpct - # remove genes that are not expressed higher than 0.1 in one of the groups take_min_pct = max_pct > arg_minpct @@ -88,7 +132,7 @@ def pseudo_seurat(adata, arg_minpct=0.1, arg_mindiffpct=-float("inf"), arg_logfc # this has the potential to be very slow. Transposeing it speeds it up a bit. # I need to undertand sparse matrices better to make it work if use_dense: - print("using dense matrix") + logging.info("using dense matrix") # extract the counts for cluster cells and calculate exp means on each row nct = adata.X.T[:, cluster_cells_ind] cluster_mean = np.apply_along_axis(exp_mean_dense, 1, nct.todense()) @@ -98,7 +142,7 @@ def pseudo_seurat(adata, arg_minpct=0.1, arg_mindiffpct=-float("inf"), arg_logfc other_mean = np.apply_along_axis(exp_mean_dense, 1, nct.todense()) diff_mean = abs(cluster_mean - other_mean) else: - print("using sparse matrix") + logging.info("using sparse matrix") cluster_mean = exp_mean_sparse(adata.X.T[:, cluster_cells_ind]) other_mean = exp_mean_sparse(adata.X.T[:, other_cells_ind]) diff_mean = abs(cluster_mean - other_mean).A1 @@ -122,7 +166,7 @@ def run_neighbors_method_choice(adata, method, n_neighbors, n_pcs, metric, use_r # useful if we are dealing with a MuData object but we want to use single rep, e.g. # calculating neighbors on a totalVI latent rep if method == "scanpy": - print("Computing neighbors using scanpy") + logging.info("Computing neighbors using scanpy") from scanpy.pp import neighbors neighbors(adata, n_pcs=n_pcs, @@ -131,7 +175,7 @@ def run_neighbors_method_choice(adata, method, n_neighbors, n_pcs, metric, use_r use_rep=use_rep) elif method == "hnsw": from scvelo.pp import neighbors - print("Computing neighbors using hnswlib (with scvelo a la pegasus!)") + logging.info("Computing neighbors using hnswlib (with scvelo a la pegasus!)") # we use the neighbors function from scvelo (thanks!) # with parameters from pegasus (for a more exact result). # code snippet from Steve Sansom, via COMBAT project diff --git a/panpipes/pipeline_clustering.py b/panpipes/pipeline_clustering.py index 3cf8c7c0..c4733e4e 100644 --- a/panpipes/pipeline_clustering.py +++ b/panpipes/pipeline_clustering.py @@ -19,7 +19,7 @@ from itertools import chain, product import glob -from panpipes.funcs.processing import is_float_try, splitall, extract_parameter_from_fname +from panpipes.funcs.processing import extract_parameter_from_fname PARAMS = P.get_parameters( ["%s/pipeline.yml" % os.path.splitext(__file__)[0], "../pipeline.yml", @@ -27,7 +27,7 @@ PARAMS['py_path'] = os.path.join(os.path.dirname(os.path.dirname(__file__)), 'python') PARAMS['r_path'] = os.path.join(os.path.dirname(os.path.dirname(__file__)), 'R') - +PARAMS['mudata_with_knn'] = 'mudata_w_neighbors.h5mu' job_kwargs={} if PARAMS['condaenv'] is not None: job_kwargs["job_condaenv"] =PARAMS['condaenv'] @@ -36,7 +36,6 @@ @originate("logs/setup_dirs.sentinel") def set_up_dirs(log_file): os.mkdir("logs") - mods = [key for key, value in PARAMS['modalities'].items() if value is True] if PARAMS["multimodal"]["run_clustering"] is True : mods.append("multimodal") @@ -48,16 +47,44 @@ def set_up_dirs(log_file): ## Single modality scripts ## ------------------------------------ +# -----------------------------------= +# neighbors +# -------------------------------------- +@follows(set_up_dirs) +@originate(PARAMS['mudata_with_knn']) +def run_neighbors(outfile): + mods = [key for key, value in PARAMS['modalities'].items() if value is True] + if any([PARAMS['neighbors'][mod]['use_existing'] is False for mod in mods]): + # this means we want to rerun neighbors for at least one assay + #we want to replace thhe scaled obj with the new neighbors + log_file="logs/run_single_mod_neighbors.log" + cmd=""" + python %(py_path)s/rerun_find_neighbors_for_clustering.py \ + --infile %(scaled_obj)s \ + --outfile %(outfile)s \ + --neighbor_dict '%(neighbors)s' \ + --n_threads %(resources_threads_high)s + """ + cmd += " > %(log_file)s" + job_kwargs["job_threads"] = PARAMS['resources_threads_high'] + P.run(cmd, **job_kwargs) + else: + P.run('ln -s %(scaled_obj)s %(outfile)s', without_cluster=True) + + + + # ------------------------------------ # UMAP # ------------------------------------ + def gen_umap_jobs(): """ Generate find neighbor jobs with all parameter combinations. """ # same infile for all jobs # define files based on jobs - infile = PARAMS['scaled_obj'] + infile = PARAMS['mudata_with_knn'] mods = [key for key, value in PARAMS['modalities'].items() if value is True] if PARAMS["multimodal"]["run_clustering"] is True : mods.append("multimodal") @@ -68,14 +95,14 @@ def gen_umap_jobs(): log_file = os.path.join("logs","_".join([mod, 'md' + str(md) +"_umap.log"])) yield [infile, output_file, mod, md, log_file] - +@follows(run_neighbors) @follows(set_up_dirs) @files(gen_umap_jobs) def calc_sm_umaps(infile, outfile, mod, mindist, log_file): prefix = os.path.split(infile)[0] cmd = """ python %(py_path)s/run_umap.py \ - --infile %(scaled_obj)s \ + --infile %(infile)s \ --outfile %(outfile)s \ --min_dist %(mindist)s """ @@ -100,7 +127,7 @@ def gen_cluster_jobs(): """ # same infile for all jobs # define files based on jobs - infile = PARAMS['scaled_obj'] + infile = PARAMS['mudata_with_knn'] mods = [key for key, value in PARAMS['modalities'].items() if value is True] if PARAMS['multimodal']['run_clustering'] is True: mods.append("multimodal") @@ -114,6 +141,7 @@ def gen_cluster_jobs(): @follows(set_up_dirs) @files(gen_cluster_jobs) +@follows(run_neighbors) def calc_cluster(infile, outfile, mod, res, alg, log_file): cmd = """python %(py_path)s/run_clustering.py --infile %(infile)s @@ -168,10 +196,11 @@ def collate_mdata(infiles,outfile): --umap_files_csv umap_file_paths.csv --output_mudata %(outfile)s """ - # if PARAMS['full_obj'] is None: - cmd += "--input_mudata %(scaled_obj)s" - # else: - # cmd += "--input_mudata %(scaled_obj)s %(full_obj)s" + if PARAMS['full_obj'] is None: + mdata_in = PARAMS['mudata_with_knn'] + cmd += "--input_mudata %(mdata_in)s" + else: + cmd += "--input_mudata %(full_obj)s" cmd += " > logs/collate_data.log" job_kwargs["job_threads"] = PARAMS['resources_threads_medium'] P.run(cmd, **job_kwargs) @@ -182,13 +211,19 @@ def collate_mdata(infiles,outfile): 'logs/plot_clusters_umaps.log') def plot_cluster_umaps(infile, log_file,): # get associated umap + mods = [key for key, value in PARAMS['modalities'].items() if value is True] + if PARAMS['multimodal']['run_clustering']: + mods.append("multimodal") + mods = ','.join(mods) cmd = """python %(py_path)s/plot_cluster_umaps.py \ --infile %(infile)s + --modalities '%(mods)s' """ cmd += " >> %(log_file)s" job_kwargs["job_threads"] = PARAMS['resources_threads_medium'] P.run(cmd, jobs_limit=1, **job_kwargs) + @transform(aggregate_clusters, regex("(.*)/all_res_clusters_list.txt.gz"), r'logs/\1_clustree.log', r'\1/figures/clustree.png', ) @@ -218,38 +253,50 @@ def cluster_analysis(fname): # ------------------------------------ # Markers # ------------------------------------ -@transform(calc_cluster, - regex("(.*)/(.*)/clusters.txt.gz"), - r"\1/\2/markers.txt", - r"\1/\2/markers", + +@subdivide(calc_cluster, regex("(.*)/(.*)/clusters.txt.gz"), + r"\1/\2/*_markers.txt", + r"\1/\2") +def gen_marker_jobs(infile, outfile, base_dir): + mods = [key for key, value in PARAMS['modalities'].items() if value is True] + for md in mods: + IOTools.touch_file(os.path.join(base_dir, md + "_markers.txt")) + +@follows(collate_mdata) +@transform(gen_marker_jobs, + regex("(.*)/(.*)/(.*)_markers.txt"), + r"logs/\1_\2_\3_markers.log", + r"\1/\2/\3_markers", r"\1", - r"logs/\1_\2_markers.log", + r"\1/\2", + r"\3" ) -def find_markers(infile, outfile, outfile_prefix, mod, log_file ): +def find_markers(infile, log_file, outfile_prefix, base_mod, cluster_dir, data_mod): """ Runs scanpy.tl.rank_gene_groups in parallel for each cluster """ - min_cells = PARAMS["markerspecs"][mod]["mincells"] + anndata_file=PARAMS['sample_prefix'] + "_clustered.h5mu" + "/" + data_mod + cluster_file = os.path.join(cluster_dir, "clusters.txt.gz") cmd = """ python %(py_path)s/run_find_markers_multi.py - --infile %(scaled_obj)s + --infile %(anndata_file)s + --cluster_file %(cluster_file)s --output_file_prefix %(outfile_prefix)s - --cluster_file %(infile)s - --mincells %(min_cells)s + --mod %(data_mod)s """ - if mod is not None and mod != "multimodal": - cmd += " --modality %(mod)s" - layers = {mod: PARAMS["markerspecs"][mod]["layer"] for mod in ['rna', 'prot', 'atac']} - cmd += " --layer '%(layers)s'" - if PARAMS['markerspecs_pseudo_seurat'] is True: - min_pct = PARAMS["markerspecs"][mod]["minpct"] - threshuse = PARAMS["markerspecs"][mod]["threshuse"] + min_cells = PARAMS["markerspecs"][data_mod]["mincells"] + cmd += " --mincells %(min_cells)s" + layer_choice = PARAMS["markerspecs"][data_mod]["layer"] + cmd += " --layer '%(layer_choice)s'" + testuse = PARAMS["markerspecs"][data_mod]["method"] + cmd += " --testuse '%(testuse)s'" + if PARAMS['markerspecs'][data_mod]['pseudo_seurat'] is True: + min_pct = PARAMS["markerspecs"][data_mod]["minpct"] + threshuse = PARAMS["markerspecs"][data_mod]["threshuse"] cmd += """ --pseudo_seurat True - --minpct %(min_pct)s - --threshuse %(threshuse)s + --minpct %(min_pct)s + --threshuse %(threshuse)s """ - else: - cmd += " --pseudo_seurat False" cmd += " > %(log_file)s " job_kwargs["job_threads"] = PARAMS['resources_threads_high'] P.run(cmd, **job_kwargs) @@ -276,19 +323,23 @@ def find_markers(infile, outfile, outfile_prefix, mod, log_file ): # --subgroup %(plotspecs_discrete_variables)s > %(log_file)s""" # P.run(cmd, job_threads=PARAMS['resources_threads_high']) - -# limit jobs because reading the h5 file simultaneouly is problematic -# @follows(heatmap_markers) +# this transforms gen_marker_jobs instead of find_markers because gen_marker_jobs creates empty marker files +# which are then filled by find_markers. @follows(collate_mdata) -@transform(find_markers, - regex(r"(.*)/alg(.*)/markers.txt"), - r"logs/\1_alg\2_dotplots.log", - r"\1/alg\2/figures/dotplot_top_markers.png", +@follows(find_markers) +@transform(gen_marker_jobs, + regex(r"(.*)/alg(.*)/(.*)_markers.txt"), + r"logs/cluster\1_alg\2_exprs_\3_dotplots.log", + r"\1/alg\2/figures/dotplot_top_markers_\3.png", r"\1/alg\2/figures", - r"\1", r"\2") -def plot_marker_dotplots(marker_file, log_file, outfile, fig_path,mod, cluster_col): + r"\1", r"\2", r"\3") +def plot_marker_dotplots(marker_file, log_file, outfile, + fig_path, cluster_mod, cluster_col, expr_mod): """ Plots some additional marker plots + Read in the h5mu file, pull the cluster col from the h5mu obs. + then pull the epxression values from the expr_mod. + """ data_obj=PARAMS['sample_prefix'] + "_clustered.h5mu" # check there is a figures directory @@ -297,17 +348,17 @@ def plot_marker_dotplots(marker_file, log_file, outfile, fig_path,mod, cluster_c cmd = """ python %(py_path)s/plot_scanpy_markers.py \ --infile %(data_obj)s \ + --modality %(expr_mod)s --marker_file %(marker_file)s \ - --figure_prefix %(fig_path)s - --group_col %(cluster_col)s + --figure_prefix %(fig_path)s \ --n %(plotspecs_top_n_markers)s """ - if mod is not None: - cmd += " --modality %(mod)s" - cmd += " --layer '%(plotspecs_layers)s'" - else: - cmd += " --layer logged_counts" + if cluster_mod != "multimodal": + cluster_col = cluster_mod + ":" + cluster_col + cmd += " --group_col %(cluster_col)s" + layer_choice = PARAMS["markerspecs"][expr_mod]["layer"] + cmd += " --layer %(layer_choice)s" cmd += " > %(log_file)s " job_kwargs["job_threads"] = PARAMS['resources_threads_medium'] P.run(cmd, **job_kwargs) @@ -374,10 +425,7 @@ def marker_analysis(fname): # #-------- - - - -@follows(cluster_analysis, marker_analysis, ) +@follows(cluster_analysis, marker_analysis ) def full(): """ All cgat pipelines should end with a full() function which updates, diff --git a/panpipes/pipeline_clustering/pipeline.yml b/panpipes/pipeline_clustering/pipeline.yml index c82b900b..a9bfc423 100644 --- a/panpipes/pipeline_clustering/pipeline.yml +++ b/panpipes/pipeline_clustering/pipeline.yml @@ -124,40 +124,56 @@ clusterspecs: # --------------------------------------- # parameters for finding marker genes # --------------------------------------- +# where pseudo_suerat is set to False +# we run https://scanpy.readthedocs.io/en/stable/generated/scanpy.tl.rank_genes_groups.html +# where pseudo_seurat is set to True, we run an a python implementation of Seurat::FindMarkers (written by CRG) is used, + markerspecs: rna: run: True layer: logged_counts + # method options: [‘logreg’, ‘t-test’, ‘wilcoxon’, ‘t-test_overestim_var’]] + method: t-test_overestim_var mincells: 10 # if a cluster contains less than n cells then do not bother doing marker analysis - # Setting pseudo_seurat to True means that rather than use scanpy.rank_gene_groups to identify gene groups, - # a python implementation of Seurat::FindMarkers (written by CRG) is used, - # We recommend that pseudo_seurat is set to False + # where pseudo_suerat is set to False + # we run https://scanpy.readthedocs.io/en/stable/generated/scanpy.tl.rank_genes_groups.html + # where pseudo_seurat is set to True, we run an a python implementation of Seurat::FindMarkers (written by CRG) is used, pseudo_seurat: False # these next two settings do not matter unless pseudo_seurat is set to True, # If applicable look at Seurat documentation for FindMarkers for details minpct: 0.1 threshuse: 0.25 + prot: run: layer: clr - mincells: 10 + mincells: 10 # if a cluster contains less than n cells then do not bother doing marker analysis + # method options: [‘logreg’, ‘t-test’, ‘wilcoxon’, ‘t-test_overestim_var’]] + method: wilcoxon pseudo_seurat: False minpct: 0.1 threshuse: 0.25 + atac: run: layer: mincells: 10 + # method options: [‘logreg’, ‘t-test’, ‘wilcoxon’, ‘t-test_overestim_var’]] + method: wilcoxon pseudo_seurat: False minpct: 0.1 threshuse: 0.25 + multimodal: mincells: 10 + # method options: [‘logreg’, ‘t-test’, ‘wilcoxon’, ‘t-test_overestim_var’]] + method: wilcoxon pseudo_seurat: False minpct: 0.1 threshuse: 0.25 + # --------------------------------------- # plot specs are used to define which metadata columns are used in the visualisations # --------------------------------------- diff --git a/panpipes/pipeline_preprocess.py b/panpipes/pipeline_preprocess.py index d85d33fb..d8c91c83 100644 --- a/panpipes/pipeline_preprocess.py +++ b/panpipes/pipeline_preprocess.py @@ -49,7 +49,7 @@ def filter_mudata(outfile): python %(py_path)s/run_filter.py --input_mudata %(unfiltered_obj)s --output_mudata %(outfile)s - --filter_dict "%(filtering)s" + --filter_dict "%(filter_dict)s" """ if PARAMS['filtering_keep_barcodes'] is not None: cmd += " --keep_barcodes %(filtering_keep_barcodes)s" diff --git a/python/batch_correct_merge.py b/python/batch_correct_merge.py index a73a6bba..b003a2a9 100644 --- a/python/batch_correct_merge.py +++ b/python/batch_correct_merge.py @@ -58,7 +58,6 @@ if multi_mod=="totalvi": totalvi_extra_obsm = {'prot': ['totalvi_protein_foreground_prob', 'totalvi_denoised_protein'], 'rna': [ 'totalvi_denoised_rna']} - # we will want to keep these so extract them first totalvi_extra_obsm_keep = {} for k, v in totalvi_extra_obsm.items(): @@ -68,7 +67,21 @@ # this overwrites the old modalities for k,v in uni_mod_paths.items(): - base_obj.mod[k] = mu.read(v) + mod_replace = mu.read(v) + if mod_replace.shape[1] == base_obj.shape[1]: + L.info('replacing the whole anndata object with batch corrected') + base_obj.mod[k] = mod_replace + else: + # in this situation the vars in the mod has been subset but the + # original object contains all the genes. + L.info("number of genes does not match base object, only integrating obsp and obsm into full mdata[%s] object" % k) + base_obj.mod[k].obsm = mod_replace.obsm + base_obj.mod[k].obsp = mod_replace.obsp + base_obj.mod[k].uns = mod_replace.uns + base_obj.mod[k].uns['scVI_gene_subset'] = mod_replace.var_names.tolist() + # merge any extra obs columns + new_obs = mod_replace.obs.iloc[:, ~ mod_replace.obs.columns.isin( base_obj.mod[k].obs.columns)] + base_obj.mod[k].obs = base_obj.mod[k].obs.merge(new_obs, left_index=True, right_index=True) if multi_mod=="totalvi": # restore the totalvi extras diff --git a/python/batch_correct_none.py b/python/batch_correct_none.py index 6acd134d..21c7b7eb 100644 --- a/python/batch_correct_none.py +++ b/python/batch_correct_none.py @@ -71,13 +71,20 @@ # run neighbours and umap without batch correction +pc_kwargs = {} +if int(args.neighbors_n_pcs) > 0: + pc_kwargs['use_rep'] = "X_pca" + pc_kwargs['n_pcs'] = int(args.neighbors_n_pcs) +else: + # If n_pcs==0 use .X if use_rep is None. + pc_kwargs['use_rep'] = None + pc_kwargs['n_pcs'] = int(args.neighbors_n_pcs) + run_neighbors_method_choice(adata, method=args.neighbors_method, n_neighbors=int(args.neighbors_k), - n_pcs=int(args.neighbors_n_pcs), metric=args.neighbors_metric, - use_rep='X_pca', - nthreads=max([threads_available, 6])) + nthreads=max([threads_available, 6]), **pc_kwargs) L.info("done n_neighbours, saving stuff") diff --git a/python/batch_correct_scvi.py b/python/batch_correct_scvi.py index f0807209..c88fde8f 100644 --- a/python/batch_correct_scvi.py +++ b/python/batch_correct_scvi.py @@ -10,8 +10,8 @@ import muon as mu import panpipes.funcs as pp -from panpipes.funcs.io import read_anndata -from panpipes.funcs.scmethods import run_neighbors_method_choice +from panpipes.funcs.io import read_anndata, read_yaml +from panpipes.funcs.scmethods import run_neighbors_method_choice, X_is_raw import sys import logging @@ -28,7 +28,7 @@ default='adata_scaled.h5ad', help='') parser.add_argument('--raw_anndata', - default='adata_raw.h5ad', + default=None, help='') parser.add_argument('--output_csv', default='batch_correction/umap_bc_scvi.csv', help='') @@ -47,12 +47,12 @@ args, opt = parser.parse_known_args() - +L.info(args) threads_available = multiprocessing.cpu_count() # load parameters -params = pp.io.read_yaml("pipeline.yml") +params = read_yaml("pipeline.yml") L.info("Running scvi script") L.info("threads available: %s", threads_available) @@ -77,7 +77,7 @@ mdata = mu.read(args.scaled_anndata) if type(mdata) is mu.MuData: - rna = mdata['rna'].copy() + rna = mdata['rna'] L.info(rna) elif type(mdata) is sc.AnnData: rna = mdata @@ -106,6 +106,8 @@ # add in test to see if raw layer exists alread if "raw_counts" in rna.layers.keys(): L.info("raw counts found in layer") +elif X_is_raw(rna): + rna.layers["raw_counts"] = rna.X.copy() else: L.info("merge in raw counts") sc_raw = read_anndata(args.raw_anndata, use_muon=True, modality="rna") diff --git a/python/batch_correct_totalvi.py b/python/batch_correct_totalvi.py index df612e48..a132b0c2 100644 --- a/python/batch_correct_totalvi.py +++ b/python/batch_correct_totalvi.py @@ -8,17 +8,15 @@ import os import gc import muon as mu -from cgatcore import pipeline as P -import panpipes.funcs as pp -from panpipes.funcs.processing import check_for_bool, intersect_obs_by_mod -from panpipes.funcs.io import read_anndata, write_anndata, read_yaml +from panpipes.funcs.io import read_anndata, read_yaml from panpipes.funcs.scmethods import run_neighbors_method_choice, X_is_raw +from panpipes.funcs.processing import intersect_obs_by_mod import sys import logging L = logging.getLogger() -L.setLevel(logging.DEBUG) +L.setLevel(logging.INFO) log_handler = logging.StreamHandler(sys.stdout) formatter = logging.Formatter('%(asctime)s: %(levelname)s - %(message)s') log_handler.setFormatter(formatter) @@ -30,7 +28,7 @@ default='adata_scaled.h5ad', help='') parser.add_argument('--raw_anndata', - default='adata_raw.h5ad', + default=None, help='') parser.add_argument('--output_csv', default='batch_correction/umap_bc_totalvi.csv', help='') @@ -60,7 +58,10 @@ # load parameters threads_available = multiprocessing.cpu_count() -params = pp.io.read_yaml("pipeline.yml") + +params = read_yaml("pipeline.yml") +params['sample_prefix'] + test_script=False @@ -78,6 +79,7 @@ # this is a copy so we can subset vars and obs without changing the original object rna = mdata['rna'].copy() +prot = mdata['prot'].copy() kwargs={} # in case of more than 1 variable, create a fake column with combined information @@ -129,6 +131,9 @@ else: raise ValueError("adt_outliers column not found in mdata['prot'].obs") +# exluding isotypes +if 'isotype' in prot.var.columns: + prot = prot[:, ~prot.var.isotype] # filter out mitochondria if params['multimodal']['totalvi']['exclude_mt_genes']: @@ -142,22 +147,25 @@ rna = rna[:, rna.var.highly_variable] mdata.update() + #make protein obsm pandas -# need to find the raw counts in mdata['prot'] -if X_is_raw(mdata['prot']): - X_array = mdata['prot'].X.copy() -elif "raw_counts" in mdata['prot'].layers.keys(): - X_array = mdata['prot'].layers['raw_counts'].copy() +# need to find the raw counts in prot +if X_is_raw(prot): + X_array = prot.X.copy() +elif "raw_counts" in prot.layers.keys(): + X_array = prot.layers['raw_counts'].copy() else: - raise AttributeError("raw counts not found in mdata['prot'], \ + raise AttributeError("raw counts not found in prot, \ store in either X or in 'raw_counts' layer") -X_df = pd.DataFrame(X_array.todense(), index=mdata['prot'].obs_names, columns=mdata['prot'].var_names) + +X_df = pd.DataFrame(X_array.todense(), index=prot.obs_names, columns=prot.var_names) if X_df.shape[0] == rna.X.shape[0]: L.info("adding protein_expression to obsm") #check the obs are in the correct order X_df = X_df.loc[rna.obs_names,:] + X_df = X_df.astype('int') rna.obsm['protein_expression'] = X_df else: L.error("dimensions do not match, cannot integrate protein") @@ -229,10 +237,12 @@ if batch_categories is not None: L.debug(batch_categories) + if type(batch_categories) is not list: + batch_categories = [batch_categories] normX, protein = vae.get_normalized_expression( n_samples=25, return_mean=True, - transform_batch=[batch_categories] + transform_batch=batch_categories ) # this is stored in obsm, because if we have filtered by hvg, it no longer fits specifications for a layer @@ -243,7 +253,7 @@ df = vae.get_protein_foreground_probability( n_samples=25, return_mean=True, - transform_batch=[batch_categories] + transform_batch=batch_categories ) mdata['prot'].obsm["totalvi_protein_foreground_prob"] = df.loc[mdata['prot'].obs_names,:] diff --git a/python/batch_correct_wnn.py b/python/batch_correct_wnn.py index ae2f373e..4a889f0c 100644 --- a/python/batch_correct_wnn.py +++ b/python/batch_correct_wnn.py @@ -92,6 +92,7 @@ for kmod in dict_graph.keys(): L.info(kmod) + pkmod=params['multimodal']['WNN']['knn'][kmod] if dict_graph[kmod]["obsm"] is not None: if dict_graph[kmod]["obsm"] not in tmp.mod[kmod].obsm.keys(): if dict_graph[kmod]["anndata"] is not None: @@ -114,8 +115,11 @@ if repuse not in tmp.mod[kmod].obsm.keys(): if "X_LSI" in tmp.mod[kmod].obsm.keys(): repuse = "X_LSI" - L.info("falling back on %s" %(repuse) ) - pkmod=params['multimodal']['WNN']['knn'][kmod] + else: + L.info("falling back on %s" %(repuse) ) + if int(pkmod['npcs']) == 0: + repuse = None + L.info("calculating neighbours") if repuse != "X_bbknn": run_neighbors_method_choice(tmp.mod[kmod], diff --git a/python/collate_mdata.py b/python/collate_mdata.py index 5b5306db..23ae8d79 100644 --- a/python/collate_mdata.py +++ b/python/collate_mdata.py @@ -34,6 +34,7 @@ L.info(args) mdata = mu.read(args.input_mudata) +L.info("loading clusters") cf = pd.read_csv(args.clusters_files_csv) if isinstance(mdata, MuData): @@ -86,4 +87,4 @@ mdata.write(args.output_mudata) mdata.obs.to_csv(re.sub(".h5mu", "_cell_metdata.tsv", args.output_mudata), sep='\t') -L.info("done") \ No newline at end of file +L.info("done") diff --git a/python/downsample.py b/python/downsample.py index a1e094bb..138d8849 100644 --- a/python/downsample.py +++ b/python/downsample.py @@ -7,6 +7,7 @@ from muon import MuData from anndata import AnnData import muon as mu +import re from panpipes.funcs.io import read_anndata, write_anndata, write_obs from panpipes.funcs.processing import downsample_mudata, intersect_obs_by_mod, check_for_bool @@ -87,10 +88,10 @@ L.info(mdata.obs.sample_id.value_counts()) # write out the observations -write_obs(mdata, output_prefix=args.sampleprefix, +write_obs(mdata, output_prefix=re.sub("\\.(.*)", "", args.output_mudata), output_suffix="_downsampled_cell_metadata.tsv") -mdata.write(args.input_mudata) +mdata.write(args.output_mudata) #This stage is the point (i.e. pre-normalisation) where the mdata file can be outputted so that we can extract raw matrix for the cellphonedb. L.info("Completed") \ No newline at end of file diff --git a/python/plot_cluster_umaps.py b/python/plot_cluster_umaps.py index 37748a15..d1a0ea8c 100644 --- a/python/plot_cluster_umaps.py +++ b/python/plot_cluster_umaps.py @@ -28,9 +28,9 @@ parser.add_argument("--infile", default="mdata.h5mu", help="file name, format: .h5mu") -parser.add_argument("--figdir", - default=None, - help="figures directory") +parser.add_argument("--modalities", + default="rna", + help="list of modalities to search for UMAPs in") args, opt = parser.parse_known_args() L.info("running with:") @@ -44,7 +44,7 @@ def main(adata,figdir): obsm_keys = [x for x in adata.obsm.keys() if re.search(pattern, x)] L.info("Umap keys founds %s" % obsm_keys) # get all possible clustersclusters - pattern=re.compile(r'^leiden|^louvain_res') + pattern=re.compile(r'^leiden|^louvain') cluster_keys = [x for x in adata.obs.columns if re.search(pattern, x)] L.info("Cluster keys founds %s" % cluster_keys) if len(obsm_keys) == 0 or len(cluster_keys) == 0: @@ -65,23 +65,24 @@ def main(adata,figdir): L.debug("load data") mdata = read(args.infile) +mods = args.modalities.split(',') # detemin initial figure directory based on object type -if args.figdir is None: - if type(mdata) is MuData : - figdir = "multimodal/figures" - L.info("plotting multimodal umaps") - else: - figdir = "figures" - L.info("plotting umaps from anndata") -else: - figdir = args.figdir + # do plotting -main(mdata, figdir) +if 'multimodal' in mods: + if os.path.exists("multimodal/figures") is False: + os.makedirs("multimodal/figures") + main(mdata, figdir="multimodal/figures") # we also need to plot per modality if type(mdata) is MuData: for mod in mdata.mod.keys(): - L.info("plotting for modality: %s" % mod) - figdir = os.path.join(mod, "figures") - main(mdata[mod], figdir) + if mod in mods: + L.info("plotting for modality: %s" % mod) + figdir = os.path.join(mod, "figures") + if os.path.exists(figdir) is False: + os.makedirs(figdir) + main(mdata[mod], figdir) + +L.info('done') \ No newline at end of file diff --git a/python/plot_custom_markers.py b/python/plot_custom_markers.py index 19174bce..0f6c24a9 100644 --- a/python/plot_custom_markers.py +++ b/python/plot_custom_markers.py @@ -130,21 +130,11 @@ def main(adata, mod, df, grouping_var, pfx, layer_choice=None): for mod in modalities: print(mod) df_sub = df[df['mod'] == mod] - # make sure every grouping var is a category - if mod is not "multimodal": - for gv in group_vars: - mdata[mod].obs[gv] = mdata.obs.loc[mdata[mod].obs_names,gv].astype('category') - mdata.update_obs() - else: - for gv in group_vars: - mdata.obs[gv] = mdata.obs.loc[mdata.obs_names,gv].astype('category') - mdata.update_obs() - # pull out the layers and set ll to None if no layer sepcificed (so the plotting funcs will use X) - try: - ll = layers[mod] - except KeyError: - ll = [None] - if len(group_vars) > 0 and ll is not None: + for gv in group_vars: + mdata[mod].obs[gv] = mdata.obs.loc[mdata[mod].obs_names,gv].astype('category') + mdata.update_obs() + ll = layers[mod] + if len(ll) > 0 and len(group_vars) > 0: for gv, layer in product(group_vars, ll): print(gv, layer) main(adata=mdata[mod], diff --git a/python/plot_custom_markers_umap.py b/python/plot_custom_markers_umap.py index e3029728..9d268f2d 100644 --- a/python/plot_custom_markers_umap.py +++ b/python/plot_custom_markers_umap.py @@ -105,19 +105,12 @@ def main(adata, mod, layer_choice, df, basis): print(mod) df_sub = df[df['mod'] == mod] mdata.update_obs() - # pull out the layers and set ll to None - # if no layer sepcificed - # (so the plotting funcs will use X) - try: - ll = layers[mod] - except KeyError: - ll = [None] - # pull out the basis e.g. X_umap + ll = layers[mod] if mod in basis_dict.keys(): bb= basis_dict[mod] else: bb = [] - if len(bb) > 0: + if len(ll) > 0 and len(bb) > 0: for basis, layer in product(bb, ll): print(basis,layer) main(adata=mdata[mod], diff --git a/python/plot_scanpy_markers.py b/python/plot_scanpy_markers.py index 74b88969..e6f167f2 100644 --- a/python/plot_scanpy_markers.py +++ b/python/plot_scanpy_markers.py @@ -59,81 +59,78 @@ # script +def calc_dendrogram(adata, group_col): + if len(adata.obs[group_col].cat.categories) > 5: + incl_dendrogram = True + if "X_pca" not in adata.obsm.keys(): + sc.pp.pca(adata) + L.info("calculating dendrogram") + sc.tl.dendrogram(adata, groupby=group_col, use_rep="X_pca") + else: + incl_dendrogram = False + return incl_dendrogram + -def main(adata, mod, layer, group_col, mf, n): +def do_plots(adata, mod, group_col, mf, n=10, layer=None): L.debug("check layers") - if layer != 'X': - adata.X = adata.layers[layer] - # run pca if not done already.( for dendorgram calc) - if("X_pca" in adata.obsm.keys()): - sc.pp.pca(adata) - # read file + # get markers for plotting + mf = mf[mf['avg_logFC'] > 0] df = mf.groupby(group_col).apply(lambda x: x.nlargest(n, ['scores'])).reset_index(drop=True) marker_list={str(k): list(v) for k,v in df.groupby(group_col)["gene"]} # add cluseter col to obs - adata.obs[group_col] = adata.obs[group_col].astype('str').astype('category') - L.debug("do plots") + # check whether a dendrogram is computed/ + incl_dendrogram = calc_dendrogram(adata, group_col) + L.info("start plotting") sc.pl.stacked_violin(adata, marker_list, groupby=group_col, - save= '_top_markers'+ mod +'.png', - figsize=(24, 5)) + layer=layer, + save = '_top_markers'+ mod +'.png', + dendrogram=incl_dendrogram, + # figsize=(24, 5) + ) sc.pl.matrixplot(adata, marker_list, groupby=group_col, save= '_top_markers'+ mod +'.png', + dendrogram=incl_dendrogram, figsize=(24, 5)) sc.pl.dotplot(adata, marker_list, groupby=group_col, save= '_top_markers'+ mod +'.png', + dendrogram=incl_dendrogram, figsize=(24, 5)) sc.pl.heatmap(adata, marker_list, groupby=group_col, save= '_top_markers'+ mod +'.png', - figsize=(24, 5)) + dendrogram=incl_dendrogram, + # figsize=(24, 5) + ) -L.debug("load data") +# read data mdata = mu.read(args.infile) - -L.debug("load marker file") -mf = pd.read_csv(args.marker_file, sep='\t' ) -mf[args.group_col] = mf['cluster'].astype('str').astype('category') -# get layer -try: - layers = read_yaml(args.layer) -except AttributeError: - layers = args.layer - -if type(mdata) is mu.MuData and args.modality is None and type(layers) is not dict: - sys.exit('if inputting a mudata object, layers must be in a dictionary') - - if type(mdata) is AnnData: adata = mdata - main(adata,mod=args.modality, layer=args.layer, n=args.n) -elif args.modality is not None: - adata = mdata[args.modality] - if type(layers) is dict: - ll = layers[args.modality] - else: - ll = layers - main(adata, - mod=args.modality, - layer = ll, - mf=mf, - group_col = args.group_col, - n=int(args.n)) + # main function only does rank_gene_groups on X, so +elif type(mdata) is mu.MuData and args.modality is not None: + adata = mdata[args.modality] else: - # we have multimodal object, and some kind of multimodal clustering output - for mod in mf['mod'].unique(): - print(mod) - mf_sub = mf[mf['mod'] == mod] - mdata[mod].obs[args.group_col] = mdata.obs.loc[mdata[mod].obs_names,args.group_col].astype('category') - mdata.update_obs() - main(mdata[mod], mod=mod, layer = layers[mod], group_col = args.group_col, mf=mf_sub, n=int(args.n)) + sys.exit('if inputting a mudata object, you need to specify a modality') + +L.info("load marker file") +mf = pd.read_csv(args.marker_file, sep='\t' ) +mf[args.group_col] = mf['cluster'].astype('category') +# get layer +adata.obs[args.group_col] = mdata.obs[args.group_col] +do_plots(adata, + mod=args.modality, + group_col=args.group_col, + mf=mf, + layer=args.layer, + n=int(args.n)) diff --git a/python/plot_umaps_batch_correct.py b/python/plot_umaps_batch_correct.py index 88dae768..5c2aca59 100644 --- a/python/plot_umaps_batch_correct.py +++ b/python/plot_umaps_batch_correct.py @@ -178,3 +178,4 @@ g.savefig(os.path.join(args.fig_dir, mod, "umap_method_" + qc + ".png"), dpi = 300) plt.clf() +L.info('done') \ No newline at end of file diff --git a/python/plot_variables_umaps.py b/python/plot_variables_umaps.py index 8c4cef73..ec44c0f0 100644 --- a/python/plot_variables_umaps.py +++ b/python/plot_variables_umaps.py @@ -8,7 +8,7 @@ from anndata import AnnData import pandas as pd import argparse -from panpipes.funcs.io import read_yaml +from panpipes.funcs.io import read_yaml, dictionary_stripper from itertools import product, chain sc.settings.autoshow = False @@ -69,10 +69,11 @@ def main(adata, mod, plot_features, basis, fig_suffix): # get bases -basis_dict = read_yaml(args.basis_dict) +basis_dict = dictionary_stripper(read_yaml(args.basis_dict)) if args.categorical_variables is not None: - cat_vars = read_yaml(args.categorical_variables) + cat_vars = dictionary_stripper(read_yaml(args.categorical_variables)) + uniq_discrete = list(set(chain(*cat_vars.values()))) # make sure they are categories mdata.obs[uniq_discrete] = mdata.obs[uniq_discrete].apply(lambda x: x. astype('category')) @@ -82,7 +83,7 @@ def main(adata, mod, plot_features, basis, fig_suffix): if args.continuous_variables is not None: try: - cont_vars = read_yaml(args.continuous_variables) + cont_vars = dictionary_stripper(read_yaml(args.continuous_variables)) except AttributeError: # this assumes that we have tried to parse a dict and nstead found a string # there is probably a better solution diff --git a/python/refmap_scvitools.py b/python/refmap_scvitools.py index 8c4dc476..d05ed303 100644 --- a/python/refmap_scvitools.py +++ b/python/refmap_scvitools.py @@ -1,7 +1,7 @@ import multiprocessing from matplotlib import transforms - +import anndata as ad import scanpy as sc import pandas as pd import matplotlib.pyplot as plt @@ -30,6 +30,8 @@ parser.add_argument('--query_data', default='adata_scaled.h5ad', help='path to query data. can be a raw 10x dataset or a preprocessed anndata/mudata') +parser.add_argument('--query_celltype',default='celltype', + help='if query has already a column with cell labeling in obs. Default to "celltype" ') parser.add_argument('--adata_reference', default='adata_.h5ad', help='path to reference ann/mudata. necessary only if you want to produce ref/ query map') @@ -40,6 +42,8 @@ help='scvi, scanvi, totalvi supported') # parser.add_argument('--generate_scanvi', default=False, # help='create scanvi reference from scvi to allow for label transfer') +parser.add_argument('--predict_rf',default=False, + help='For totalvi and scvi, if reference model has random forest classifier, classify cells in query') parser.add_argument('--impute_proteins', default=False, help='if totalvi reference is used, allow to impute protein levels') parser.add_argument('--transform_batch', default=None, @@ -56,7 +60,10 @@ args, opt = parser.parse_known_args() sc.settings.figdir = "figures/" +sc.set_figure_params(figsize=(8, 6), dpi=300) L.info("running with args:") +args.predict_rf = check_for_bool(args.predict_rf) +args.impute_proteins = check_for_bool(args.impute_proteins) L.info(args) threads_available = multiprocessing.cpu_count() @@ -65,6 +72,9 @@ params = read_yaml("pipeline.yml") query_data = os.path.basename(args.query_data) +reference_architecture = str(args.reference_architecture) + + mdata = mu.read(args.query_data) if type(mdata) is mu.MuData: @@ -72,6 +82,20 @@ sys.exit("we only support querying using RNA but your mdata doesn't contain rna") else: adata_query = mdata["rna"].copy() + if "prot" in mdata.mod.keys() and reference_architecture=="totalvi": + if X_is_raw(mdata['prot']): + X_array = mdata['prot'].X.copy() + elif "raw_counts" in mdata['prot'].layers.keys(): + X_array = mdata['prot'].layers['raw_counts'].copy() + + X_df = pd.DataFrame(X_array.todense(), index=mdata['prot'].obs_names, columns=mdata['prot'].var_names) + if X_df.shape[0] == adata_query.X.shape[0]: + L.info("adding protein_expression to obsm") + #check the obs are in the correct order + X_df = X_df.loc[adata_query.obs_names,:] + adata_query.obsm['protein_expression'] = X_df + else: + L.error("dimensions do not match, cannot create the raw obsm counts in query") else: adata_query = mdata.copy() @@ -87,7 +111,7 @@ raise AttributeError("raw counts not found in query, \ store in either X or in 'counts/raw_counts' layer") -reference_architecture = str(args.reference_architecture) + if reference_architecture in ['totalvi','scvi','scanvi']: L.info("running using %s" % reference_architecture) else: @@ -136,30 +160,132 @@ vae_q.train(max_epochs= max_epochs , plan_kwargs=train_kwargs) adata_query.obsm["X_scANVI"] = vae_q.get_latent_representation() adata_query.obs["predictions"] = vae_q.predict() - df = adata_query.obs.groupby(["celltype", "predictions"]).size().unstack(fill_value=0) - norm_df = df / df.sum(axis=0) - - plt.figure(figsize=(8, 8)) - _ = plt.pcolor(norm_df) - _ = plt.xticks(np.arange(0.5, len(df.columns), 1), df.columns, rotation=90) - _ = plt.yticks(np.arange(0.5, len(df.index), 1), df.index) - plt.xlabel("Predicted") - plt.ylabel("Observed") - file_name = "SCANVI_predicted_vs_observed_labels_query_data" - plt.savefig(os.path.join("figures/", file_name + ".png")) - + if args.query_celltype is not None: + L.info("Query has celltypes in column %i, i will plot what predictions look like from scanvi model" % args.query_celltype) + df = adata_query.obs.groupby([str(args.query_celltype), "predictions"]).size().unstack(fill_value=0) + norm_df = df / df.sum(axis=0) + + plt.figure(figsize=(8, 8)) + _ = plt.pcolor(norm_df) + _ = plt.xticks(np.arange(0.5, len(df.columns), 1), df.columns, rotation=90) + _ = plt.yticks(np.arange(0.5, len(df.index), 1), df.index) + plt.xlabel("Predicted") + plt.ylabel("Observed") + file_name = "SCANVI_predicted_vs_observed_labels_query_data" + plt.savefig(os.path.join("figures/", file_name + ".png")) + if reference_architecture=="totalvi": + # temporary fix is disabled for now, need to modify to allow to fix query to match reference + fix_query = False + if fix_query: + L.info(" will do some manipulation of query data to make it match to referece structure") + if args.adata_reference is not None: + reference_data = os.path.basename(args.adata_reference) + mdata=mu.read(args.adata_reference) + if type(mdata) is MuData: + if "rna" not in mdata.mod.keys(): + sys.exit("we only support querying using RNA but your mdata doesn't contain rna") + else: + adata_ref = mdata["rna"].copy() + adata_ref.obsm = mdata.obsm.copy() + adata_ref.obsp = mdata.obsp.copy() + else: + adata_ref = mdata.copy() + del mdata + # match query and reference vars + + #adata_query.layers["counts"] = adata_query.X.copy()#already taken care of + # are the following necessary? + sc.pp.normalize_total(adata_query, target_sum=1e4) + sc.pp.log1p(adata_query) + adata_query.raw = adata_query + # subset to reference vars + adata_query = adata_query[:, adata_ref.var_names].copy() + + adata_query.obs["celltype.l2"] = "Unknown" + adata_query.obs["orig.ident"] = adata_query.obs["set"] + #query.obsm["X_umap"] = query.obs[["UMAP1", "UMAP2"]] + + + # obsm slot slot for protein counts can be called either prot_counts or protein_expression + # in the reference anndata and model. allow a + + if "protein_expression" in adata_ref.obsm.keys(): + pname = "protein_expression" + + elif "protein_counts" in adata_ref.obsm.keys(): + pname = "protein_counts" + + if pname not in adata_query.obsm.keys(): + X_df = pd.DataFrame(0, index=adata_ref.obs_names, columns=adata_ref.obsm[pname].columns) + + adata_query.obsm["protein_counts"] = query.obsm["pro_exp"].copy() + + + for p in adata_ref.obsm[pname].columns: + if p not in adata_query.obsm[pname].columns: + adata_query.obsm[pname][p] = 0.0 + # ensure columns are in same order + adata_query.obsm[pname] = adata_query.obsm[pname].loc[:, adata_ref.obsm[pname].columns] + + else: + + # For now make this only work when reference anndata is available + + # if args.reference_var is not None: + # reference_var = pd.read(args.reference_var, sep="\t") + # else: + # sys.exit("I need var names from reference to make sure query has overlap") + # if args.reference_prot is not None: + # reference_prot = pd.read(args.reference_prot, sep="\t") + # else: + # sys.exit("I need var names from reference to make sure query has overlap") + + # if args.reference_prot_assay is not None: + # L.info("name of obsm slot for reference is: %i" % args.reference_prot_assay) + # pname = args.reference_prot_assay + sys.exit("need a reference dataset to check for matching query entries") + scvi.model.TOTALVI.prepare_query_anndata(adata_query, reference_path) vae_q = scvi.model.TOTALVI.load_query_data( adata_query, reference_path, freeze_expression=True) - latent_choice= "X_totalVI" + latent_choice= "X_totalvi" vae_q.train(max_epochs= max_epochs , plan_kwargs=train_kwargs) - adata_query.obsm["X_totalVI"] = vae_q.get_latent_representation() + adata_query.obsm["X_totalvi"] = vae_q.get_latent_representation() + #remove this after finishing + adata_query.write(os.path.join( "query_temp_check_totvi.h5ad")) + + if args.predict_rf : + L.info("predicting celltypes") + predictions = ( + vae_q.latent_space_classifer_.predict( + adata_query.obsm["X_totalvi"] + ) + ) + adata_query.obs["predictions"] = predictions + #remove this after finishing + adata_query.write(os.path.join( "query_temp_check_predictions_totvi.h5ad")) + + if args.query_celltype is not None: + L.info("Query has celltypes in column %s, i will plot what predictions look like from totalvi model" % args.query_celltype) + df = adata_query.obs.groupby([str(args.query_celltype), "predictions"]).size().unstack(fill_value=0) + norm_df = df / df.sum(axis=0) + + plt.figure(figsize=(8, 8)) + _ = plt.pcolor(norm_df) + _ = plt.xticks(np.arange(0.5, len(df.columns), 1), df.columns, rotation=90) + _ = plt.yticks(np.arange(0.5, len(df.index), 1), df.index) + plt.xlabel("Predicted") + plt.ylabel("Observed") + file_name = "totalvi_predicted_vs_observed_labels_query_data" + plt.savefig(os.path.join("figures/", file_name + ".png")) + L.info("Mapped Query to Reference, now some utilities") + if args.adata_reference is not None: reference_data = os.path.basename(args.adata_reference) mdata=mu.read(args.adata_reference) @@ -175,8 +301,10 @@ adata_ref.obs.loc[:, 'is_reference'] = 'Reference' adata_query.obs.loc[:, 'is_reference'] = 'Query' - - adata_full = adata_query.concatenate(adata_ref, batch_key="batch_key") + #expect the batch to be always encoded as `batch` in both Q and R + adata_full = ad.concat( [adata_ref,adata_query]) + # add param to decide if to recompute the total embedding or to recalc the embedding + #if "X_totalvi" not in adata_ref.obsm.keys(): adata_full.obsm[latent_choice] = vae_q.get_latent_representation(adata_full) if args.impute_proteins: @@ -200,10 +328,10 @@ adata_full = adata_query.copy() if int(args.neighbors_n_pcs) > adata_full.obsm[latent_choice].shape[1]: - L.warn("N PCs is larger than %i dimensions, reducing n PCs to" % adata_full.obsm[latent_choice].shape[1] -1) + L.warn("N PCs is larger than %i dimensions, reducing n PCs to " % adata_full.obsm[latent_choice].shape[1]) -n_pcs= min(int(args.neighbors_n_pcs), adata_full.obsm[latent_choice].shape[1]-1) +n_pcs= min(int(args.neighbors_n_pcs), adata_full.obsm[latent_choice].shape[1]) run_neighbors_method_choice(adata_full, @@ -220,24 +348,23 @@ model_name = ''.join(e for e in model_name if e.isalnum()) -file_name= "umap" + model_name + "_" + latent_choice +file_name= "umap_" + model_name + "_" + latent_choice fig = sc.pl.embedding(adata_full, basis = "umap",color=["is_reference"], show=False, return_fig=True) fig.savefig(os.path.join("figures/", file_name + ".png")) -umap = pd.DataFrame(mdata.obsm['X_umap'], mdata.obs.index) +umap = pd.DataFrame(adata_full.obsm['X_umap'], adata_full.obs.index) umap.to_csv(os.path.join("refmap/", file_name + ".csv") ) file_name= "query_to_reference_" + model_name + "_" + latent_choice + ".h5mu" #change this to mudata mdata_save = MuData({"rna":adata_full}) -write_anndata(mdata_save, use_muon=True, modality="rna", fname=os.path.join( "refmap/" , file_name)) - - +mdata_save.write(os.path.join( "refmap/" , file_name)) +L.info('done') diff --git a/python/run_find_markers_multi.py b/python/run_find_markers_multi.py index 66225d17..387e40ab 100644 --- a/python/run_find_markers_multi.py +++ b/python/run_find_markers_multi.py @@ -3,14 +3,11 @@ from anndata import AnnData import argparse import pandas as pd -import re import numpy as np from scipy.sparse import issparse -from itertools import compress - from panpipes.funcs.processing import check_for_bool -from panpipes.funcs.scmethods import pseudo_seurat -from panpipes.funcs.io import read_yaml +from panpipes.funcs.scmethods import find_all_markers_pseudo_seurat + import sys import logging @@ -47,20 +44,14 @@ help="testing limited to genes with this (log scale) difference in mean expression level.") parser.add_argument("--mincells", type=str, default='3', help="minimum number of cells required (applies to cluster and to the other cells)") -parser.add_argument("--use_dense", type=str, default="False") + args, opt = parser.parse_known_args() -# args = argparse.Namespace(infile='prot/ncomps3_nneigh30_meteuclidean_dir/neighbors.h5ad', modality='prot', layer='rna', output_file_prefix='prot/ncomps3_nneigh30_meteuclidean_dir/algleiden_res1/markers', cluster_file='prot/ncomps3_nneigh30_meteuclidean_dir/algleiden_res1/clusters.txt.gz', testuse='wilcoxon', pseudo_seurat='False', minpct='0.1', mindiffpct='-inf', threshuse='0.25', mincells='10', use_dense='False') # script ------------------------------- -def main(adata, args, mod="", marker_columns=[]): - # check the X slot actually contains data - if adata.X.shape[1] == 0: - L.info("skipping %s because adata.X contains no values" % mod) - return None - +def check_log1p_dict(adata): # temp code to get around this bug # https://github.com/scverse/muon/issues/40 if 'log1p' in adata.uns.keys(): @@ -68,97 +59,135 @@ def main(adata, args, mod="", marker_columns=[]): L.debug( "log1p dict is Empty, refilling dict") adata.uns['log1p'] = {'base': None} + +def load_and_merge_clusters(adata, cluster_file): # load clusterss - df = pd.read_csv(args.cluster_file, sep="\t", index_col=0) + df = pd.read_csv(cluster_file, sep="\t", index_col=0) adata_shape = adata.shape - # add clusters to adata.obs adata.obs = pd.merge(adata.obs, df, how="left", left_index=True, right_index=True) # check check we have not lost cells in merge if adata.shape != adata_shape: L.info("some cells lost in merge, not all cells have a cluster?") L.debug(adata.shape) - if issparse(adata.X): - bool_list = (adata.X.sum(axis=0) > 0).tolist()[0] - else: - bool_list = (adata.X.sum(axis=0) > 0).tolist() + +def get_header(): + # write out the correct header before we start, need to werite it out before, because we use append mode for the multimodal find markers + # and we don't wantto end up with a new header half + if sc.__version__ < "1.7.0": + markers_columns = ['cluster', 'scores', 'gene', 'avg_logFC', 'pvals', 'p.adj.bonferroni'] + elif sc.__version__ >= "1.7.0": + markers_columns = ['cluster', 'gene', 'scores', 'avg_logFC', 'pvals', 'p.adj.bonferroni'] + return markers_columns + + +def filter_zero_count_genes(adata, layer=None): + if layer is None: + if issparse(adata.X): + bool_list = (adata.X.sum(axis=0) > 0).tolist()[0] + else: + bool_list = (adata.X.sum(axis=0) > 0).tolist() + else: + if issparse(adata.layers[layer]): + bool_list = (adata.layers[layer].sum(axis=0) > 0).tolist()[0] + else: + bool_list = (adata.layers[layer].sum(axis=0) > 0).tolist() adata = adata[:, bool_list] L.debug(adata.shape) - run_pseudo_seurat = check_for_bool(args.pseudo_seurat) - use_dense = check_for_bool(args.use_dense) - - clust_vals = set(adata.obs["clusters"]) - markers_dict = {} - filter_stats_dict = {} - - for cv in clust_vals: - L.info(cv) - # 1. set idents - adata.obs["idents"] = ["1" if cc == cv else "0" for cc in adata.obs["clusters"]] - adata.obs["idents"] = adata.obs["idents"].astype("category") - # check we have enough cells - if args.mincells == "None": - min_cells = 3 - else: - min_cells = int(args.mincells) - n_cluster = sum([cell == '1' for cell in adata.obs['idents']]) - n_other = sum([cell == '0' for cell in adata.obs['idents']]) - if n_cluster <= int(args.mincells) | n_other <= int(args.mincells): - sys.exit("not enough cells in cluster, min cells: %i" % int(args.mincells)) - L.info("Cluster %s number of cells: %i" % (cv, n_cluster)) - L.info("Other number of cells: %i\n" % n_other) - # adata = adata.raw.to_adata().copy() - if run_pseudo_seurat is True: - filter_stats = pseudo_seurat(adata, use_dense=use_dense) - L.info("number of genes remaining after filtering: %i\n" % filter_stats['background'].sum()) - adata_rg = adata[:, filter_stats['background'].tolist()].copy() - - sc.tl.rank_genes_groups(adata_rg, layer=None, - groups=["1"], groupby="idents", reference="0", - method="wilcoxon", n_genes=float("inf"), corr_method="bonferroni") - - markers = sc.get.rank_genes_groups_df(adata_rg, group="1",) - # remove adata from mem - adata_rg = None - else: - sc.tl.rank_genes_groups(adata, groups=["1"], groupby="idents", reference="0", - method="t-test_overestim_var", - n_genes=float("inf"), - corr_method="bonferroni", layer=None) - markers = sc.get.rank_genes_groups_df(adata, group="1") - # filter positive only and adjusted pval < 0.05 - # markers = markers[(markers.logfoldchanges > 0) & (markers.pvals_adj < 0.05)] - - markers.head() + + +def run_clustering(adata, + cluster_col=None, + cluster_groups='all', + layer=None, + methoduse=None, + pseudo_seurat=False): + if type(cluster_groups) is list: + # make sure the clusters are strings else rank_gene_groups crashes + cluster_groups = [str(x) for x in cluster_groups] + adata.obs[cluster_col] = adata.obs[cluster_col].astype('str').astype('category') + if pseudo_seurat: + L.info('running rank gene groups with pseudo-seurat method') + markers, filter_stats = find_all_markers_pseudo_seurat(adata, + groups=cluster_groups, + groupby=cluster_col, + method=methoduse, + n_genes=float("inf"), + corr_method="bonferroni", + layer=layer, + arg_minpct=0.1, + arg_mindiffpct=-float("inf"), + arg_logfcdiff=0.25) + markers = markers.reset_index().drop(columns='level_1').rename(columns={'level_0':'gene'}) + else: + # find markers + L.info('running rank gene groups with standard scanpy implementation') + sc.tl.rank_genes_groups(adata, + groups=cluster_groups, + groupby=cluster_col, + method=methoduse, + n_genes=float("inf"), + corr_method="bonferroni", + layer=layer) + markers = sc.get.rank_genes_groups_df(adata, group=None) L.info("Rank gene groups complete") - markers_dict[cv] = markers - if pseudo_seurat is True: - filter_stats_dict[cv] = filter_stats - - - all_markers = pd.concat(markers_dict, axis=0, keys=markers_dict.keys()) - all_markers.columns = marker_columns[1:6] - all_markers = all_markers.reset_index().rename(columns={'level_0':'cluster'}).drop(columns="level_1") + filter_stats = None + return markers, filter_stats + + +def main(adata, + cluster_file, + mod, + mincells=None, + layer=None, + testuse=None, + pseudo_seurat=False, + output_file_prefix="markers") : + # check the X slot actually contains data + if adata.X.shape[1] == 0: + L.info("exiting because adata.X contains no values") + sys.exit("exiting because adata.X contains no values") + # prevent log1p base error. + check_log1p_dict(adata) + # load clusters + load_and_merge_clusters(adata, cluster_file) + # filter out genes with zero counts + filter_zero_count_genes(adata, layer=layer) + # filter out clusters with less than min cells + if mincells is not None: + min_cell_df = adata.obs['clusters'].value_counts() + if any(min_cell_df < mincells): + exclude = min_cell_df[min_cell_df < mincells].index + L.info('exluding clusters for haveing too few cells: %s' % str(list(exclude))) + clust_vals = list(set(min_cell_df[min_cell_df >= mincells].index)) + else: + clust_vals = list(set(adata.obs["clusters"])) + all_markers, all_filter_stats = run_clustering(adata, + cluster_groups=clust_vals, + cluster_col="clusters", + layer=layer, + methoduse=testuse, + pseudo_seurat=pseudo_seurat) + # make sure all files have consitent headers + all_markers.columns = get_header() all_markers['mod'] = mod - all_markers.to_csv(args.output_file_prefix + ".txt", index=False, sep='\t', header=None, mode="a") - + all_markers.to_csv(output_file_prefix + ".txt", index=False, sep='\t') + # subset to significant markers signif_markers = all_markers[(all_markers["p.adj.bonferroni"] < 0.05) &( all_markers.avg_logFC > 0)] if signif_markers.shape[0] !=0: - signif_markers.to_csv(args.output_file_prefix + mod + "_signif.txt", index=False, sep='\t') - L.debug("test_debug") - excel_file_top = args.output_file_prefix + mod + "_signif.xlsx" + signif_markers.to_csv(output_file_prefix + mod + "_signif.txt", index=False, sep='\t') + # write out to excel + excel_file_top = output_file_prefix + mod + "_signif.xlsx" with pd.ExcelWriter(excel_file_top) as writer: for xx in signif_markers.cluster.unique(): - L.debug(xx) df_sub = signif_markers[signif_markers.cluster == xx] if df_sub.shape[0] !=0: df_sub.to_excel(writer, sheet_name="cluster" + str(xx), index=False) - + # when relevent write out filter stats if pseudo_seurat is True: - all_filter_stats = pd.concat( filter_stats_dict, axis=0, keys= filter_stats_dict.keys(), ignore_index=True) all_filter_stats = all_filter_stats.reset_index().rename(columns={'index':'cluster'}) - all_filter_stats.to_csv(args.output_file_prefix + "_filter_stats_" + mod + ".txt", index=False, sep="\t") + all_filter_stats.to_csv(output_file_prefix + "_filter_stats_" + mod + ".txt", index=False, sep="\t") @@ -166,48 +195,28 @@ def main(adata, args, mod="", marker_columns=[]): L.info(args) L.info('\n') -# get layer -try: - layers = read_yaml(args.layer) -except AttributeError: - layers = args.layer - # read data mdata = read(args.infile) -if type(mdata) is MuData and args.modality is None and type(layers) is not dict: - sys.exit('if inputting a mudata object, layers must be in a dictionary') - -# write out the correct header before we start, need to werite it out before, because we use append mode for the multimodal find markers -# and we don't wantto end up with a new header half -if sc.__version__ < "1.7.0": - markers_columns = ['cluster', 'scores', 'gene', 'avg_logFC', 'pvals', 'p.adj.bonferroni', 'mod'] -elif sc.__version__ >= "1.7.0": - markers_columns = ['cluster', 'gene', 'scores', 'avg_logFC', 'pvals', 'p.adj.bonferroni', 'mod'] -pd.DataFrame(columns=markers_columns).to_csv(args.output_file_prefix + ".txt", index=False, sep='\t', mode="w") - if type(mdata) is AnnData: adata = mdata - if layers != "X" : - adata.X = adata.layers[layers[args.modality]] - main(adata, args,marker_columns=markers_columns, mod=args.modality) -elif args.modality is not None: - adata = mdata[args.modality] - if type(layers) is dict: - ll = layers[args.modality] - else: - ll = layers - adata.X = adata.layers[ll] - main(adata, args, marker_columns=markers_columns, mod=args.modality) + # main function only does rank_gene_groups on X, so +elif type(mdata) is MuData and args.modality is not None: + adata = mdata[args.modality] else: - # we have multimodal object, and some kind of multimodal clustering output - for mod in mdata.mod.keys(): - if type(layers) is dict: - ll = layers[args.modality] - else: - ll = layers - mdata[mod].X = mdata[mod].layers[ll] - main(mdata[mod], args, mod=mod, marker_columns=markers_columns) + sys.exit('if inputting a mudata object, you need to specify a modality') + + + +main(adata, + mod=args.modality, + cluster_file=args.cluster_file, + mincells=int(args.mincells), + layer=args.layer, + testuse=args.testuse, + pseudo_seurat=check_for_bool(args.pseudo_seurat), + output_file_prefix=args.output_file_prefix) +L.info("done") diff --git a/python/run_preprocess_rna.py b/python/run_preprocess_rna.py index 1df9811b..84283e3d 100644 --- a/python/run_preprocess_rna.py +++ b/python/run_preprocess_rna.py @@ -234,4 +234,6 @@ # save the scaled adata object mdata.write(args.output_scaled_mudata) -L.debug(adata.uns['log1p']) \ No newline at end of file +L.debug(adata.uns['log1p']) + +L.info("done") \ No newline at end of file diff --git a/requirements_minimal.txt b/requirements_minimal.txt index 47337598..40a54d2b 100644 --- a/requirements_minimal.txt +++ b/requirements_minimal.txt @@ -26,6 +26,6 @@ scrublet scikit-misc mofapy2 pysam -#mudata @ git+https://github.com/scverse/mudata -#muon @ git+https://github.com/scverse/muon -#scirpy @ git+https://github.com/scverse/scirpy +mudata @ git+https://github.com/scverse/mudata +muon @ git+https://github.com/scverse/muon +scirpy @ git+https://github.com/scverse/scirpy From b4fe99ea96be2f1c7bfcf587a376e3bb3016b823 Mon Sep 17 00:00:00 2001 From: Charlotte Rich-Griffin <25774602+crichgriffin@users.noreply.github.com> Date: Thu, 9 Mar 2023 11:24:00 +0000 Subject: [PATCH 2/2] delete cache --- .../funcs/__pycache__/__init__.cpython-39.pyc | Bin 305 -> 0 bytes panpipes/funcs/__pycache__/io.cpython-39.pyc | Bin 15439 -> 0 bytes .../funcs/__pycache__/plotting.cpython-39.pyc | Bin 8006 -> 0 bytes .../funcs/__pycache__/processing.cpython-39.pyc | Bin 11442 -> 0 bytes .../funcs/__pycache__/scmethods.cpython-39.pyc | Bin 9646 -> 0 bytes 5 files changed, 0 insertions(+), 0 deletions(-) delete mode 100644 panpipes/funcs/__pycache__/__init__.cpython-39.pyc delete mode 100644 panpipes/funcs/__pycache__/io.cpython-39.pyc delete mode 100644 panpipes/funcs/__pycache__/plotting.cpython-39.pyc delete mode 100644 panpipes/funcs/__pycache__/processing.cpython-39.pyc delete mode 100644 panpipes/funcs/__pycache__/scmethods.cpython-39.pyc diff --git a/panpipes/funcs/__pycache__/__init__.cpython-39.pyc b/panpipes/funcs/__pycache__/__init__.cpython-39.pyc deleted file mode 100644 index 8882564aea2c9024831cb90b1db6fe1942e71330..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 305 zcmYk0y=ucS5XU9?BOf#m;-wlALWWW@cj%JMi*RMdQIM@eClk^)>MMEe)K@43dP;}( z4*l=<;}BA_S}hsX%k8EO&CewNOOWj8;D(}QEqAQvHRp*1v8b~XNuP zu&VsxbhYA>u%2=W!E1wPhmLH4x;OFUTO^r-^D)xmwVF*Vq^CnhF$sW>Co6uK6UW{{6HzE2iG;q4g5luk(j_2)eD=}yZ%>Wh*Dj46K{d0L72aQrGJJ_hy VBV4EXQ}5f+ndkHiEGZs%>L1~|O<({3 diff --git a/panpipes/funcs/__pycache__/io.cpython-39.pyc b/panpipes/funcs/__pycache__/io.cpython-39.pyc deleted file mode 100644 index 757b8d63cb336e2331653b42e4313dcc1881eecd..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 15439 zcmcgzTaX;rS?=5P^v=%C&aPG~tt87W*@`sIW>2vwd-{h>6K77Nk>zXDs zp|`ZQ?&&(^jh5jV{4_n2pO$CwGw0d-%zJr$I-Y~4*_vt>yh3~0o96yjtJt3LX4)mM z)SmTb+hwoZp7Z9~^WJ=h|A4o^{SSHv`MKyF!ZX)8++OmQ0Aq`M>qz^kceH)XJH~yS z)_v{c-tnye{hkZ>sg~6~;hjL;lVV1cUe>*nqR`Ytp<#+yQAYUz@qWMHKe(hJJ%w~y z%r9x;z-y*h5C>n*dG8QUp>7d%hhEdgVX=g|)8e8yChkM;Gvc_o|7Fdq2v?lI^C9sx zpdSGAgMdCI-T~-$iZ%a)e_Wh?%@Aiq1vSg!nPJUC;+?2@m-v8K5v#zsBF>73foWAd zBHoSXS@Ebihv&oMG4UQeAHl31Z@s(ysP||#({tWA@sr}c;tBMAO#GNQ52*Ku9~T$! ze4N+rz21ApPl)%S{Dd$*tW}=PF2@P-~_6OU=lawN}MW zi_))&>ceMmRNK*wbb46YI+spyVe4E?WZB5S5v2~-*28OAZc8-qnp}OW)(pN4GScv- zxq2t4*P?1o)S_Co)eNJw*bMN>(67@F&>{b15L`qnvd-mNv~{Vu%cAzZdWgayje5yHN-WCv zEiv8C1(qoGC-C+);LkN^muA-G9OmO4l(p;ChF^<%(ywmRtL<7On>V~8jb;mr=8JpV z(p%()Ce|gj+kUjw0S#s|YYMMW>vsJ>RKrN7&S2fvg$v(nv75B9pCV!t|?0#>DU0B;IUBGGimuyot7bw{`Kd1JvS7)SKO`Vo|E zw3KxD0H7ua`XDN3-|P5>UpF-npU@K!tk9RVo5n@$Rr94k?OTztXC`38OCZN3keqx3 zWS+y<4!(eT)KR_2fk+!5eS1lhjmR4HNOX}G&TaE!Mn8u>jj=vc=<~&_&lm6NQvlp_ zhP!y%%6gw6yu#YVXts0oUR95y_Vvi#%M8Y`U)94eCpz%_^>D6v%rzf4gIS0J*{ls|W7Z3pHFB8s z6mT9G%$n*(<2ldJeb_F|Ow_To6-=8i4CL^v}xuY2*)QxQdKx z8n}vLDZ^#TCxGh+a8WL!HF0#u0-jGtg}rHUjF=j*EM-_`kiRd>m&Ebmn7)yX*+@A%@GyAjR{^6UAtjGPs9;P(VBP;jI+@mbqVUI_%9`0C=N6{lI_vg6BT(E#Q z&EQ+&H=6x9Jj>$55D!tpiEIw@+-F{#L|sPTVA2P$ucrJ=azH#lrNBY%L9QP{^%p>g z2a<)vPV$M9OeKY6Iw>YINhz64%E??Z-zYz?>Dtw@Fu~MM#TN@MKF|08WNZYM>b5$N zu`CQjce}FzHYr`dR^M`4&A@k?fgAQV+RZS8Xk3YIMBH7HEKQwF|3;N~<$7AA?10(H zMQSx=P-6*FG30%ykjIgDvoUN$E-g}V*f<4=U(fJURu(a+_=OicJ-1#9+)mKiRcjF; zkCg++Do>%)_E|=u3hNm)+y07srQ_B%H>JOcMP7D8-yf4w+_I9+@cV|ZFEwOsv+V~_ z$U{mGll0uE$$CflRfrkctY<9?UB9+dtpyNjHPToZ-YXAbBf+zX?UmK2-Cd2Rst{EW zdF@`D<42{K=2*6-xp1r2^;4rOQmYwa!ydq#k!^MAsUdxczurcghfL|UgD}nEjet2q zJ9d|kovYYop~8p2DUO4DdXlQpTD-1X+UiMFfSkexP-BsGFSK1euDQF4ZohbVa` zCGVnSg_5(BJWR>EDS4EVbCf)W#G4^Pf?lZRn&zQk&{8ttu&h-HiNW;XrXN%lod*G{ zeIWFZUJ`x*iDoV8^AM0WMB_1h%X;ChqPd8glCfy0mLqympEr*{hz?qddKn_rG9E%p ziE{cVw@5fpU({W08>J+6vp<9T0iUfG3; zQ1S>jNx%DsndTZnt?j1{w1{fE*9lSw`%r5&qusg*u)&A&4!(pn{3en}+k;$BG^I4g z<;d8E(wVRd`4KHRW@%Aw_X!>R14<}KQmBx+_M#RTXEjK0v!}^Rk&|dQ^$#mtR$?ez z;O8X!ThY`W#jQDdX#Q>sSfp^&OtQe;1$)lLm-3#k30taf@>TSVpLn>Q} zSJ4&{9cfm{FhHrGWGIA=QNp*No=^=P^)=;lfCO)zBtNj*pG5tH>PU7GeVM!($Y9HL zkjO%!?SWiqi-L7nI0iOxInt?4n1>+N_hE0Kp6sKGS$%G-p2+3}Rj*(nfK6o?@s>O5 zYwBD)i+8HS@l%we5y$WO5qGI#rMc^pOwr4PYsB+V&(cypL-ik^gjVA2C7KXB-%WGPK=?N* zW@>^%O88~8)Q&n;nADXV>^ho@Hv?YV?LibvzX3~dPz-wwDmsu=HKJP232WDUVoZ0h z&6HXjyy0)~f{<8}1}<0x&B>1=Unliq@G+4YQWY&XnM8=<*C8a>Sen>gpM$8+%>4z{f`|etj2mo=mrzbfXR(;a*5Q zB-+AZZV5W$ENodnqDp1)jcG^BVwRy?sSg|*wGdd?LL-!5&h9{Cf?gFyIIS*t1RJQe zKyV{19NW5M0E0RlLei5)vNK+l{?%R+M!2g)G&r8ix*7&f;`Z**K^40I3~Ax{TC3+X zBLUL)p~{2n@Aj?hQLw-}37_)sZI1D#zcWUr0)lrZ4D@2^IpOCu=Cz4XVvP zw1Gnx$#+m7siR~J=a?~o6_DHH69Wk(K)tXk9|r`iN@5*X|7+N0aT}Z_V{Igza?E00 zCFh)Q#1#GtV)`{3WPsDK=j?0Obom9A7sW(h1bd>oy(uxXukRaRPS2CA0g{2GON>Tn zU|!0qDydo}+n@&z<)oO!sQ-g;mIrb)DDlYV24yC;`E^!{HLBn<1Eq1yGMU?S%17Ru0IaW z#7Kw4`FTEawlh7kZv>aX3rH1caV>lY`#NeLLoI2PpW@nLR_mgc)Pv7+ExZB@Tj2A^ zuG55EnNIvDo2}R3YlBtg)+FgsE_5?KB4NELySr|)asje=J6s9fFh5fnn^UyYb>{E(RV0mA&J*%8_A_fpx2-4@4Q8J=KqbXCP+3`lcbqw26zIBB^E2vj^%F#}50jaXe z;{u^)qSxINs*tFq`5ypR{7)?E~zNr6-4rewT>+tn~8Yw5tTbXSGRn+R~@;i9u zz}3G48DX{Htz@YDfGQ(vF9#2dvkzY#w@eWxc;xcrl@Y}yQa1jRsB2Us3grLXb0Wa@ z9qu)gIGGjXB)RQ^D4bWLDJ9x(Yrn%b6vNFXIR(e1GQnp3@_WI zb#)Ugx3{rU@3dFL`Zc;$fs>>YtZuYA8<1$wdZem=`q+t-ADgi$cM^G7*6y5Cucf0m zIt8v$J3|`Qsc`)YJE7u2z0=(#z*AwW%Xn^J*721NPF9#YgDkxDe0KqAQ`Z?`tZse? z4E^JkCk3s-Ck7op0iiIV_YL}~DaU1{+<}v@i3tp|N@2eSB@%K96$ir|$Qv3T?x%1q z0Q5DvIicW^vP=5X0iteVuL*W|2Yl-deVvck8+xkWKvK+Bzq;f5U3K`_z2Tg1)!G}P zcHuNo(ZSQ=6Jfnu9ox#4^YG|kmKTswu<6X;PKRVp7#~(kPUFig=CB$2m_lV%+CZro z3|Vy$U0c1V*&-zYcM#SZjxKyUupMer)^<}Hn+%k&CvXE|VkvNQ*6M_QYVl>8*{*d{ zi*~ItC9e|ZJ|T4m>yuhopfI&}uG3RqA&6f<;@Mm0Fh)O}8f;G{A`^n(g_yl-)PW>nCZ!Ib2$0)9g5k*PIgKc7G$V?V|i`!SjQ zSbRu_Oh+A5`7sz^>X@35#m6)!Y#dW-J16pE$JE9#wGxX3)-S*P!Gz-Y|9diJh&ivx z2mCGEtB5ap+- z_{)^sqU1A3;3pf>l>8N{{VXL+MLd8JwUqC9FOf?L;6yIu05SBFZ#ky?9BQtRi#zbE z;cj~aMpB)Q#Nb2XX&PTLbi%(v@{0b7(SYwm7%w{!Sbt=^f;!Z|_kkLi;&|e6j8mIa zxYkhPyjJMjw+um7V3$MzmpJn%kt@s-Qz+4`ZJvabU;~DP0UZBquG}QyLoGw$bIVAa zZF@Tpk@Xd{IHJVxsN0*`*C3XrlYCPbv$yrV>01y=iBq)`2h>-N`%HbL{m=zKE!ZUi z0OkoTy4&qYV)`4qtUIw+b*MJ6=XXH4FgO6GY-#YGmg>YE*pHwb!t6nZfxs4}(!1`V zFRG|gH>hHAOqVbwJ|y?Niaa78#+U7HdKE7Rga(#FX`F4eLWw5@i(jKqPvZa3 zkA|O>S;8-W^ZQ}H{G(xKE%^sy|LGqMJE`CD2V;M3fW1;+yIpQjzd9v?l7R_)?t7{C z^49=D`V1IGQz$oDJWbT36csl<8f`&8z(1Mj1Az;}nT=w&rb9o92p06b!X7#7BzbNs zMRxZ}VmXqCRTe0@L#?`+@nMtY;F6&U7e@_tB(&d+mJgoVti#Q49qvG=YweNy(9J%|$x5`~X5I{# zZ(@QU7};_#>sdfy6U*^Zu@B~x=!djlHLik9j*dV0ps{jR{wx5YRNCqVBgT^!&Xb-7 z1-^UX+Taw#=dwtX@mUTdf0*Nr^T$E|Ndp;sxM;~u~_|m}}GtF-wtxEGd@dJm~CgNicym3rlgx?+aYj4@`h#b-7mw{(O z<0IlCX8shuVT43EI!Jy-Mr48zteg5?4mpi-`yjo(5!q3G8)jevvrCZtH*f1vVH*xb zxZdHsxJ_hN+@9GhCFr*ghvTrs?njhn_sR)7I};Ka_Ig}_^PFymTwG%_E#bx*A{j!J zoM>=M49IGisSsl2af}uZj~+&l9IXcZG3+@i$8J%p})1P+D5TQ*n)F$?)8z%+3x z06V%jpdF@rgFOvSVBFf!&A}iap`L=9kiuBLi56kv0!4$l4M!C16X&suln*{WhUW&y z-?6)4h~uul!OGaXzQpi0>tQOKX!NcOzaW(+AO9k2QBiFN?1E$j8<=#ju1rtnt8c4C zEi5D19BFczX7aHLg)co0Di32Asomg#9AEPnrJ_vF)9!S(oP7cnF#Q%>R_O+5G9c4(yYU{ zqld6yt8Vm^e~z{ZdLRpq6iGhQzlZ*GrfFyQ*?HO2aBHL~_WsgeBUp=t1R^kTb)cLw zfWY-O-(GGAlWe;4(BTyZ(3B@^)aB1>2ylXJM>#-~CjV(J!wCe98xfay&SbquaOB$j z@(`{x+0NW&mbnj&+2S0I5iyKJxj)GWl zxfVoVZbbLQcHzl)*^L~AwO#T`(JkIbpLzO)W%tR8SC-u?Pd&HnzMqmOuRQtGOGAu& z!+^64cT*3iau{eu_@D|R6bo{baFLtC9bOAIn!Za8l_4H~Wpl-KfA%Eac=8GNqXbRB z%PzUi`1bJROUpxOnhQQBF|-VrQBOekl|zr<)-n`yV=}PyXr#EgTt&PUP@EjOB(QOfxRH*IvsNU{1-~)B{%rLBqloKC!SFl4vS7E7O1R1{NiJ_>GjF$WD zUI9KbISabEqEV6<F%X@4k0Q247sKriR67phAkD6W zsIkWQ0rpz{yHx8+ zv%}jQ7e=-PM?%E2&yWK8^bM>ZNA4-b{1LFfE3BxmsW6EU74q|xkj@^LpA4erMn~+r zPr<*2r@ivYw$r<7N=mGIV{qaT&^QhptulWS$D$_!w(tbZmyA zpVoxPiM+K{SPj2NFMop)QuScoj7%Al7l=oQ{C21=HpmX}4h~IxB^jW%TQjK*$*>HTGL2+^NGgiBJ*wm%s@@ znQyl_c4ERDB0G=_LY8b_#tPjzt(9?OwA6v~=YkvWM; z{C+RLui-ZeW_F!nMg{IMO}u>HKrl)1Ch4s+pOHQ9X1eT*n)qJwj4r zZ&t~yk03pzUQVx%fOEgK8T3(RY7ZysJ#$cya* z!!g--H~snt{^yuGk}kfj|0b+;Gx7UBOoM8N^@e-&@*@T@%x**t>!(D6sxF#AK*X&*6f&-^gUSEp~>RN�uX$$ z5zmauh!9AP>l^WG7OE3cED$Q?3q%j5d=ZZb+l9=$>q;9@h0w)BQ-%~5+D;42kW_RZztNDtgl?&c zlKd;ezC=a+N@{iRTZq)6OZT+EluJ>?@a{{$-NEh99}=7cEg!+J(EUF~eQ#01E{rcw z?zaiV60*~+ILc@^j2bV;Lqu$(l%CbVIJ{hger4i0gQb!;sQZaw!WKqZg&xswium^< z;kW3>J_I2ASc2oFvC z3QC)lSD+oj;(+$_NTyNZ(lRn`8uB67j&MLAXp4UF1Ncp?VIe#P_FtyKl@q;*Yan@< z!bix}z+oiDx+EhOhSgu-L0dq9Lbafd>lM;i(gP#Cg?@rU!449JwM9%hbpro-m2>(| zs$aWnapXS(rEXw77Lm~q+e$2q*bl7<{e;BQ8D?VcP|Zem7HylF^PbN0UdFu1Hbx{? z1U`W{Ea=dr^<6E3OHs&AQJ_bHkhQDG;Y1kgskII@3l_s?fQ{!V|A}(nK$0#%e76v; zz>^C$s{uem#H;B-{tE$7T*Te;`Zv_g)4U4EP_UAUn2=OIF&ElL^5ioo{U$QxSh}f4 z6va57z(NTAl(`It$l`|^bXlx}`>pqlxWw1O@T2#2gwC!Y^ugYCb{jglEL)-ND@r38XIj^il*9qko!jeka<=pN2~ zj_a%#gQt8RK=LC5sk1W{6tq^+*DKZ)W`U~kV_-SZzdLP-7i;oMl$@vJ5E2?Y!Zbs6 zbf{Q~PjNLC2-=Gz5|>09o*m++0Nq`65k)6nxntQA%RX76FwdeJZxuI6>J0TFWls4` zrYNQ~>fjPd5~-yHcmnIBrQt<0gpvGYUOsaQvquJj_q{gyzf8?^2f|_c92_oRK{2&9 zfR6)pIhvG%C^-;M1^e(F0N*!!jyjVbr4%M|nkcP_^a})R)|=g3XqAeIe~ns+*YTm} z&6X(@Xh21AsD{W6$3HcET9u}eH1%!A(NJwaPfaNMJ1LBDXz$+8J|5W#Q z#vwh^r_SZCQ|Fxe|Nr~`zrAv~sNwU|PoA`2e@WB+nHr}*1{$~UP5uc5*EkEc22)qP zVW_Lww;DFf$DD?XK0Pe-i;be{Gs2mEsZpAYmm6huPo+`8)eL9*)kbwPHrJTrC2l zD{1l$YHrt0o9(_gK;0N7e#3aX+kNc3c&+p1o2}QS*0$-M$;ZXDgzrz_n|u#Nre|yo zC1*EOxu%V^BW=u%w2__JV?EQ3uw@#@h>_V`Khn=go7^~Jr+a!<&(3hnU7efUIy{l? zIh=8uJBRv6&+u}Fn)q7k_6l@}bbN=#ujnIZWTeI33|6K0O1wlX!?nDISDTrY&P6M7~^Q-aNn#&UL@VNm-^jU?7fSla$>=K z;qyCnUFxkpS@4p!AMqgS*7al3D=~vYvK`5`9~si=1l=Uq_hq5uwSCMP3t8}aH?Il? z*<&j9ofo+wulY6Nr*{)=J(%2gM(iwQ6pQgTKiL|_Mu;s}@7)1U~ zz|$>hhkh69k-E1ljkwdP&q{k31@~~HoyK91q|z3CKfaHNgg*$eDrtL3x;OBpfsM*a zLVA(!IPeoWlf)uz-uJ>G8pa?@rJD>}bVO(-5eYy#^#yfF2U{HW(QA6UK_Uw@SsJ9F zUw37xE8^jx+1`qSw%;g1a#F{Z=)`$P8s2UqZ91l*ke0f!QI69sUra>tOq}Q=O}vg{ zdVF`IJLn{@ZS44AxY71R8h3>^*xK8GK8s{yKk2>k<8N#vqP;QjqCqh5ljdO93WKNZ zjm|J?CmZzq6z}lLVDDr<%{cOB@sy;7LSrSp%8L3Tn`h_site&S_8haB`MEqlgW&yT569Sx}&QfF<#s@#cyOz&*d|` zbjU_dWM@tn*NHKG%*$OyGmi{j8M$cztHINIFsE^eR&-nY;31w;7!@(1#&{?;qjEYsuBLNY@qqE!6%DoW0b5ayA)U|4FxHvO$y^w#krmcx z{257JJ%kU=TwZ)r8@czsfA0^u1{0mzY}i9MuI9DL_0?lKcx8Bkyp$+YK`)^Djo+rv zpWlSF+FQv$Y5R`4ig}Dnb_X_XDVr7+-9^Uum}PZcT%+czDC!H+<{}kOYHkJydUG<9CAZferr&}?XR?ADibT8v&k(_i67d@D z5Lf9oBjE5My+%2|eQJiev*VptB@4tAk|xW_iyx!=megaK6C-b?J29Mrco{><=TK;@ z%3PR6mA%9&Fpi3DGZ%(pvhy&J68hKKB_@6Wl5KL2(+|Bd`hE#tH4!{A421k}%6h<$ zj6uoa&jXej-)3J$KR595&HdZAg||~X<#cL0LAq5N#6gtS$aRRM))8^P=GDmOrarIT zh6}wTyuMH46Xj_^NpYR^8ZBnS|hw7wd$|sO6WLl5>0AQ1R12+y0=)2a_iANmLmPy~O&Ep1S z;TQb~Is{dAMVxwR5Jz%hD;B|i9Hn09g+Vv!Lj)Ra?}lE>53#jb^1!=k(B4)(iJT?e zKiy$KqRzwxJmOQE%Qd>W+azy(8xn~kR&aPa)8YEWucEI_kTv~Koen^4F+v3n51Z_< z`>@F9Rca9@nDo3t4sk^rv1nC|7$+mQ;ONwdI--d?92-Y;Jm3dw+O{rwsWrBN>1;V? z8$@-SO=*b8+c!~Mm*xN`wQo^*U77a&`pZczNNQdf);7J(OSQHaMR8gq0H`}a%vK_P z6swkIkdVD6!@jhVEpOn9cc|O)g^1}DSbf6fh67t5s_^W|t--$h_BIVQNg117wAb{K z#M@H}Pre6`f>W$QQL4JDpNI8`0KIKekm-j`yMV7c%?{KjrSmiobb<-6ER*boRObO8 z1a8#q7`F*TEH$cF4Z%lh;V$by$5o>%_>jR(l;BxCoFLd08`Oqajb9Z z3OWc_3cX@xgHgG>45~WAXAcQ{oZbvZO0)u8%S2u=Hh*aF1-^J_L8YBur8f&8`FP}J zZm*iT<2k7O5?>~`LeD&UK#rgGj2DlzkC-T7Yzbp%jX`+;(d}fWFZAkCY@OHr4HkJK)n}~aU@IAfK?<(x?}AX9bdgIW z1JVxKoPe>$NuYeairADuegU=D4)w8>Ax>SLd%P|3dF=$ zXp1E<#<(r*JbV(j&?g;)Y6Pd^cc_JM1`!qFx15|B{U8!=(?B_C$HH$WULRJBu*csW zh%aCOft0?xE53!Q_+1iYzbuVy#H7Sm8uMU}O`f0UG-z`s?PuR{v5DK7?V-5u|0eED zh&pTRB1qgNyf47f@1w1J!5K^~C571%`RLh_}6)zL1X?I0|2xO)!hruvEsiwap0 zEU^Ien&E6zN@p^gf{~Bd2d^Of048z5t1>UHXnZC^X5i3*lv=^4g0+@2XHlD!&I;53 z>NA^`vso3=^bqz`M{@{q<}#=SFRNU`{AeL7^=7kqUda}4cNKF}c4027u4+I%!aFys z;whMg);n5+eD$<^R?QYC@(Uz??|YD#>&*i-7m4;|rSSq^%4gs%LMB4VTmmRKkyylO zPGpW0Z^m8A*(@;vk{z!cX+YfPMoZaJZ*{!JSJepB*ZCUNKVtX3dGBl4%=jE3AA+%x zQtkEqFBfb0ya=+){aVXQz@X!(ww0!XSS8>;6a56@ryLYITX^)mxXnxzCS`AT$Mu<`r6M$$DiOsiJ*e=?V;aEkBOy8yM+`& z+^?U{0mnczlH0r4_JdHorDOnIJnl;~@qI1~bSYcULk}B?jiKxzeIeHl3^m9cU=zL#P@iTDF- zp;1scM_W*DBlllFOPl){Dt-+Ga%4fbnWwReRLyBd9>Y?2B3ab1qb5z2q&h+3ERl8R>+(Kg`(R{?UX zh;Qdv4G?ybKu#y1BL`Vx^GF1Vzr+1!AdB>qV2fUG5~BD}`%s6UBA)I()O)Z?cm-qs z^$%Hx5fgvl0+J8}Qa{Bu1Xz?-;Wnb2JMcWJFEigmIB06Av4>EA{s43irIz0Zgx%Gy zsnT_*1OP}TKlIE6cqaHUmzl*?D!O)J&`*@vFWTRZdA*xp`s}EDm*Yl z-WLHe^1kLnwb-0$*V()IS`D=_Wz~}1hh^I8iIYg37AN{B{){T0y{Z$It*EqeCi`XV zvqsTjhbVX565q!h(h5?)pX3Hz{3#-oKS?p6M3f1ehtxmPlZcRseNMWTK=D4Ael#@v1Rb4bts#Qd?S_O7`>lTI0v2NEE2MXE?7`T?*9u= zC5==mWl*I{P%}M6m`E^v0A>W}ajAu|TQrs-W2pM5ZcgiWHPOSQbrbePd15n9BYT4Z zWtuG*V1j*v>LB@Lfa-we`vcXFp$JFRWzK(!HANPOH%fO5>X*XVwA2(q*4ub zrIM+xLUwCp@0alF1+oky`#mjMg{_%;S#E2Yjijmlz82|MwbV-O6e-or_$n+`afPoy zQyrB)MN0Dkcul`fjPd(~T8&HVJfeI?0UX7JWtVW|f1MS0c?CMnD{Gn>EA`5~3b80! z8A^&{)hLpnvwRj@#=#z{hc571)yK10LA|Drz))<3-*ca~e=@AQ_$C)oD8pJGF&r9} zo-!B~k}naopgW0Ta}CwSzd>63I~ChhoDF1^QT;9TE2DaW{OL3emYg#8r-ntI`Ham| zWOlMNiZRvFO4Fs~cAhVWDSP|Ob^QHw9Z#)Sdw6>h8c_D4Icz z0K8xeVVapBZLC7GdR3Xi(`2)X0A%IuG!;Q>nEJ}|CxP0>HQ(%<52t(`c(`%4Ij%VkV@~YF_WKEQ#Hz`8a{8h1kXHx8+?ssQQ)ix zyLk*!z>h?e5cjOg$nqH_`5WAtTtcA{mf}a1tGW&^K!5AN*Hv7RTZYr4KX`~A+D_FK zPat16HuG#f0b8BnhfVy8Z%;D~b13`74dxJxa$f&$?>&Xm4Xsyg_Ks~4AxtFD8;`QqggE4syNZpoebKiQ89ZvX%Q diff --git a/panpipes/funcs/__pycache__/processing.cpython-39.pyc b/panpipes/funcs/__pycache__/processing.cpython-39.pyc deleted file mode 100644 index cda0969931630816365f844b8edf02fe8b1aa752..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 11442 zcmb_i+m9R9d7l{$$su>OUai-QPzqb%c|`(Tc?tfnz7By@(g#$ zC5P0Rp;pUsbFq|K0je}ivo=b6m5Z`&qdLP`~?L%uLTPDr7s2gU?6ROzweOT zrB(@ypxm7^=YIVzzwbK}95_(4@O$&;Kd*lK^Op7ROzi!+NW6?E{;6$Q%2IaIYT0Gm zX548yWk=qbaz@^6*_C&;oRxR3oRfFH?8&=OF3NkVJdJmzS!f+7ACSIW#AnKnNj%$} zZ5=EhY|WMDB%NzM-kLAZOFZ9vqIIZzNa9{|s&%-0SmK4|k=D`jQHd9u$68O8pOpAi z^LXp2@>7UUs{?nf@(DGg9>e=-HLDKd?W>o=!{IY(?juJ%uIBG%%O{b3Hhfk+@li${ zQit!ljykqxt0$FJbJ6PfN49!Ooj|M8>S^VpzccC?brSEh>RELP z?~-~>oyPll^<`C3&!dmes8`iFwSd?;wWyX*wxF)63+e^rEUMSkMfD%`u-Le$zKC&`)k~P$i)+>zPNJ4kSJao#?z76dZk4Wm%rh+6Y4Pni z)YlVTkJeU@%vIY_HAoQ4wnDuYmYme9wp$$?#&Mdht_Ssqwc4w3Mb!{1#NE|6Xmy$) zv%lKC8YBT?`B$Q77|rb#>ri#8VHI>T{xR<1c@(9btOOV=&eO6;wZ zw#C19u#W52dC&z`;PzfA`uN@K)~g$S62^(2wEZ{?{d&z0{8*OcSA)o34gI96qfq(v zC}{_Nt=SF|X?x?^3Kp1^jI{XmpxF&y)4Hv>7Sq84OV46fpS-=a)~Ur8mo~#@bEz8W zq`jtt&id9;7uy#v_2R~*&tF=Kb#dTYu2`%SPgsN^sOCx+vz*X zK8h{eaI}B3eA7kR8Kz(E+k=d&cl|Dz(EXjFMP72p88Dee6!MzLd0=s_S7GhyF`f z{BNFN-kEd$nWW9ol$?3bmsa&C?$kjy{_2*$o+O?4^3qbA1l5i9tx#jpQsXn_=T?ALU=<=5)XFhY*sjZ~=p zWIgn2?Iz9%1n*bty4r2w_#&{tW$EXQmABr$;n&#YvJWA1e(_SJ+U&;YaWT1_EUxwZ zw)U6D!`vc|Ua2exQPA9q!9B@VC%m%EVg0yUtJQB`+0*a1S$x?)iQZ32*>Ajk1A`(5 zgPELHyIO-#LC+v&)eht6OyY00^#)}L*B~8XPBJ`6d$sM?V?PRSi`0>2?oMIRKP%L3 z)-zm3Nyh%fem|XLXa}d8!*iWSeFVW8nD@O|GBu7IY`Nq7*tzLq`|Z@OKwLRG)Q3>- z0iuiNQ~UPAL|0sH*5iZ(euYxvzwlVS^O!16tTNt{3J(n?Y+;1y`=2k|s0s7C-HTOhaMI z(v*G*ZPN_2jxbH?a=x1TmmZMy)Nx9|YBNl8pb5m&ybe3fpc;xO)d!?Gq(CVTxlgjw zr|7f?!F)Kbu(soQF+xXK+L{Rp$h5}}X8J*xrTltAt)dk5;WZ!4#Zy+wXz z;OWc!c=DzjnfOG!?-ckB=BB%gSV!XCLGvg3}k?XcZd>%Ak8BgMwlV4CTE2ATtM z0-ECM_pHHz9VlUe5Bc7I3prnh_74*uaDcoaBG^A`6n;#`twLSKC1`~OIRi>@>^-DK zST!0r2{tGhCORP405|&WhzuF|*oy#KrXI#VPh0djGvjj{0>X7~4KOT=$ek;Y?~5|c zx5iRRk_+RXBxTDl`v$y74siv4CJVD`jOP2vt%*DivkGk#szv{e#BbHt))V3l9Rik6 zZ9!Y=a6yLwqL3q|Cx5dJRrxK`@TF^DyzpDB*$OrQRZwlbM^w0!LyIdPlNfqufoGz2 z;_Gg-Td9KG!)D<>9PPx`1oyHj=|{MiqC-T?%(e|KhVE}zM!*aiIV}*~ku_p{0+qqq zFW|My&H5Q6LuLshpB>MZ6RE6jRd4{0k|?E_)T_s4)v31`=m8&!ikQ`Ktn^j*Bokca3gegBR}xarqCyP%a8}bfkoj!DuZnBS0i2g zD=;CMR@G~@5NfX)3eA5SP9$v#w|w*Vdwvips(WxVIoM=xy+=A0#)Tn~z*-oRJs2|b zZ_7Y1H~`BLgtx_TSoGg%cm2%(9FMl9M`;B)JErGxxAh3{K#`MTAyET6bQGNPSG$S7 zg|1ZFCy%g4v$Aoos5b6sG(U92xp8S2P#XOIiHq9?<>(hd0KYOz>Pu)XtY`>VNSOLY z81dVT2qT`xsq#?wQ1YYLa0OcQ3v6_m!50~j&-9lVyv*Pg247}C{9yRzDqbyJHl1>>}Bgi>qEN+O{ko^VH`J4sGhge%WMgzGU1*U>}< zcHdq@mq+&sS}6MrXfJi^Fa}GGei}{7S;4+QcG*C1KpwVurP^*jBz%8}y5`JJsrF`M z6mCVDJnF=WD?|jQ(;8e0)CUZnOd8@2y5u=r#}n_>5irN?25`h49f8@FI0GUYLzhAV z5{&8c4_V0lC>8Y^c%>PvB|W0U7EHHFH-h1;D#!}g+8U}Il6oB1l$V+?Q;?m0o-vM? zy47}PtCZKj&V&<1`m2%)t6!x#SwGYVcCMV?2*XZgwQ4r*E9`(XOJ^V6ukRz*pfX~y z%%RK)yo>g1W}b%GZ?R1VV&ew1AQeiLM?Ly6&WBtQbIyHN*`(FJHL!rd;JWBTcGLbQ zjsm@9L2qu~YB(l$*-C)1fYDsgGq6KhnBzs-QerppLQ zSd1`kU|lLwa`+_OEcE=lmu&@iTaumr1~gxLzW z#JfmOZ~+a3EFRekl23xgCk~cJg@$dQ`FY*a&!PxoM97Er(<~@*E!3Pf44lGhyWOk^ zrhjN=h0XPPb-hxL)IXtWOmeksT|qbjlSgDM^O!*?=Yfr%7=K8<93Q0oK3=H@Jpaqa zMe)C&-e}wa@rm_w+zcKW_gjc4am(~E`W;7qxes3!%<1FG*|u-Ki&O@w=atJ;;ta4$ zgz$;0vVH3icI*Q?Iz~TOl11ycNDL^;jmpF^lV3^Q8@e0bNV9? z&``(5ov>Q3ZTSle&&7Ta`)onfEKPcE?6&YfT@(K-IthFplYsA`pRH8Nr9}!dW8A|) zOJI35+iPpU1`u&kt(_KbV0z+SOI@}vIpt~4VmE;_B7}fA>?Vd^$C4zSfjg+WQK_|c zg>8O@ofRY50xm0xlMIeBr)bY0MaPQXK#hGTGiFjkz~1c0!>H_S`x{nt+ISZeyD>Zo zoD9w>^04P%&JXN&tjIYJF$k7)SM+)Gk-ESsVee2R&_Z_9xCUx%I6lCXv@n;V z+bCuOLA8*^F5Q*7!QzL=e2yg*)~#U4(&xXSno*uM-Gat3m~&=q-9`yKu!5p0b#-k* zk&?Q{9Qjiugq}c+@mvB$DaIhylw=xk4aKLC2M+xJAQa*$+sNIr^_3*AGJVI;BgL2G z-asy#io;yi_1HpXmqD*c6<~P>b!KA|5un&jOuN5TK zb>l!>4U^3IEHIu^?*tf& zt{}7pj2sUGzhtx9isIhkiE>pN6_kssfv&byh?9dGBdjGv1n^tJ=|~@Ye%QfEX<92B zpQdM6uOoo*qte1v?8-InifC8eWr6pVaLaI57J@Rot*OydSLt>qbqSY2nMQT0S~&-^ zh|uDjtdonvO?dPlCdYZX1gOwnHSmJi@;}fZ_7PZFvLMVNQ3^$S8bmoK%;Vc4 z`B4h3|CRL`Q^=tvh9~vDOJb$1HJBdt(KtZg>3sk<{WGbb1%f287JwbmzLHI529Gf| znB5Tp9?c{N2Xo2eDo+H->m`u0x6*qGUVo5>p^|+Rrf*^cv_??MKxl#zD(xusG@S{E z=@=z9O1AFd-Ft~$E#3+=y9?k6sK*23d+ajH!s|Kel#8zXotO{4V5c%tH zdF)yk=@hS-zF|y$=n9vs&~#(z?=tuvg0ja<$UL1E0SVWJy*?x$rvjx|r_xFkyirzs z3QaB11+>qn8FS77tH5smpIr91$V)ic9sRqgw+}0fbv7{zjuJlC9Sc_}I6fS-9H7B0 z8GQ#e^zR^m$Q`X2@PYhPK0x-7UD!qIyGm?!=$zoRrP6Uudy+%+Lmx6!Ugy?K$xV*wqA z0}(<{1gd~gq4g;<8g4G&&7g=ym9KNc(s_Yh$wfn&mvP`2 z{zFJOc#Ke+cR#jOw;;>t%U%6}Tq4{PSp@qSc7}MG452_}OCMR1jSUML;^90S%Oh!?n`U+7{CO2RF6R(&!|~2qokO79joTN02U}+ z5B>rICHq{*GeMvreb}<}8%)}>8Do%$gYYFZ8k*y&-FiqNipS;RJ}MRG*CVNPALDRd z!Z;@p*>x3=btrmnSsES+GSQ^HJ)WBI5a59m50&+pc`lEpwQo!^tLy~xSpS&82?oMTW4S4ty*wI2 z_+9@2s`ZMO$C!6wt*9#8BWopGq8-lL5zvXV!mk$S$H(U7f|!<%KE**N6^uTEyf;UK z_MRT{^!^hbUQ0cR>iZ@^xh<36eV)QbC%C*FU^*aqAwEv&fwS_7?Zb!sN2sJfU_f#4 z3n%rJ(WD;OsTL} zD*h=ezQBsJxTC3zk8A`p!DQ0cN31s?K+R^K$9m&)D5Q3uf!Dq008}Ex-5JE`NT5z< zy_^r=dCED2>va9ss6V0QU1ac~4=8>fLEqMgr~_+u9Zs3Qw0{WU8Q~&rGN+GkOpN9? zeaiaaxp;scQV-6nN~;T$mfOJX#ajGGTkhaqeIR@Cn_QBI9%#PK;o?vK4dN5zLsY=& zzKTaIT~3vE|8lEOG7@lt$qi5g)U!)Cq$6zVdA!mrZjDs@5N#8jscz>T&VplD*hu}i zcu#B?_0Q-L1tyw7JgY|E@Ea3wRA#>)l$}0~M=di*=2h063@#NOZC0(kkD1AI+5 zy6-F~ujMRWN5nvs2mgQw^gn%^^CPef$OyZMh(q!( zi_;H&&@~(WOkjg~alLJ4qvuT;A)G3(!OEs4QgA$=yqIJMIee935Dd`AAio2o^4^Je zo$dU#N87Fsy$_qsqu0Kl>*ukb9e%3Bt=@%-+K3xi4?_DwW zb;vKl;(hx~3#t|dC=79{ydK=5$rOh$(ARN^%PZz)uq7Y4z;_l0b=-UnYfHB)z0b(q zZH)Hg>t|YA5PG4H>&!LW_rhynI>D(BWTaV18ei^TAmfp=Cu(G`7FEz`G*5^Y3`3Rd zPwWMBO8&=pN^quihMWbT`{?dq_voKNweOAD>O`wsf$d1U{!vT?=fg7Un+dzBUWHBA zt9Oik5p}aBqvErRR;&qMHD!mP)SURQVh}UcK)Seg}gr@lc_K8`w;3mcGn0$*3m+OxQ|b@7?i5fJ#it^XP2B^#gotq7RF zCUuZ=4e#&ZEdj*9?@2B5x|G5+fO?Zp-(pSS48wRUHkPfdpP0&gm6DCRtuCe9Zw>VQ$rAbI4!id`F3MB`xril?2YsSX0=Tk6DT_Zu&c0i#QtD)pB8&(dn3r z5i0XEkFNpfi4xyXnlqO>dWuB_er{-boo26d4309mz@Wt78Uq3VxrQ@$2=6jUx4I4( z(48uNM)Bo{4XtUhiaHc!D=JL1h{#dnOQ0B$?cQ@{V+u^YWgK(28D7DuSuss1`j966 zV94|PIEMIRD&tLgPk6`Xo}4>ceB#(?&&A(C??CZb@wn#}XNxDjyqEQ|v)=y!rDmr3 diff --git a/panpipes/funcs/__pycache__/scmethods.cpython-39.pyc b/panpipes/funcs/__pycache__/scmethods.cpython-39.pyc deleted file mode 100644 index e952e5a5c46b14dfc2b7f4335c52dd4914734c9a..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 9646 zcmcIqS&SUVdG6|KdX7D~T$0PfQd5>mjzum>>t-mCvScc9q>T{~4$+p<>YeJD?b)90 zVRaApaL2oU5!AQB(}g1iKhmG7^b znLQ*K$U}NCHFf;;SJhwt-4jer)igX`{ph3CGf!#SzfxoGW1(>!-|QCCG_G+LYaWxo zx@X{5kIkOtS&YVv*zQ$4N7l`_+N*grS!Z$GtIIhJuYq4Hp6X3|)4dsQW;8zQ&C2mZ z-W-1I_)zb#cUX>9;`!bY@5l)5sCN|OPJFDl;4Pr9%4_lQ-U;u7tk>g{y;I&PS#QLj z>pkK<(tFf<6zwTKeOL2bKEr46dzv2#ALDadtgi9He7n)WsyeF{pU20wBPx7ZQ?>uL3YRwlO63v&ItaPF*>jxqWQMc3O%;#;? zW|CC&f;id?{XETsxT%-UulMsPO@g?rydLGDKqqFrF?c1&1Jo)nCCRv&1$V-B)Om;- z)G$;P5$}ZkIL*@O7Y50=o}~E&40S^8avVs`&1vqPP-Hmv7_L1?T0{nSEfuSo%OerC@_55V4`ElN zbbCQAqIGw9kh?&U#p&9j#3+dOn^vhO{nEg4r4_FCdzVX-wo<39d>W7)%A1uewZu=tb zHoKU$dSSkj@@%ocQO;mRpO*90UXKwmbrhN{7EpeI;uUmCil0M;D2}g=uTg;3oGqX% zw64k7aqT3;m2N=N|F8ak=`X_D0}wxR1DA(M7LI6s-px|C7j;&0cNu+QFBKc^${^{8 zjcam4U#2@z*D~mb!fgjxzV~j#42~&gQItkD=#>WORa)6f0N{BVqPtNAk*=Jc&{J&# zeTtkFaSHX$4}U<L!dH^{XJGSO14%(#$g{qxZ?U1Z zr^k>zHX!uE;wCrW)!RBCtyQyYaEsgbS#B2A7CWg8D-z27OekAIIk~lM=ap@zs1#0H z-(lO;9bDtET2#9=ULlOb8tU~?y{n=RJ> z#O#57i52~OuyhX(=Zd+)z^)JT8n55g^7-O$FT1ALOA6IcWiP&>`6b1|GHBp$L;!=(a`Lb4hORNF_^(F?VyEb}?P< zUJtS?*vKCLjir~aPgL;J6TAYXgSv*aLFlDMj>9|QXqUTds1rT2=#mX`vpg7`@!824 zmxI=7N5J%OcR9$Q@1eSJ7{A{l@}q!tQ0JGFULic#<7i>18L=7)bsNtujxRwjM1Z%rz#f3s!%8gvB?Sf;w>S*AlF@^mfhgIJ`tG&+vnm|LE`b~VZM^JR2 z)!sh(mN8^QeP|RcT_|*Xjo;I6n{Vzh4o$Pv7cZ48DA{t!S_iF0?P?rlIf=}RJz)JW ze6`KfO89-H3l^`&X)B1c7f~6*ZaxZ@qEn2SO0rx#%1B6w^P08NB&wHGz&`uOQML9(&~#S$;Dy9FdOvy@rED7 zY0~lKRgz6DD|srfu{46^OnNP){?g^r*o^vR^-duC1WRF0r6pO?s|BLt12@j*O^@ny zAZkxh8?IQpMLk{>&Xq5bdxyrjfQ9^XO* z9Hw>QAm7mh9KfNGe7Rr+>*^Gfbj=;O%|lZ{Sl}kcb$pEb=<&%V#q8tSHe6yEx1eM&#|NUm`p_IdH5w1Eg9bt03MR=NfiCpXi#V&sbt(YhI-Blotwwu z0*iN4t*fh`qH7Exk=s^bAnv2|OXXbyVFxtAyX4wKj%gl^ybBFLj(pc1TFLeOSpKg# zh8FIuQULSN&Yfg7&3PCea+c_d`BBftCbH z6xbhO-HeQ(26f(m3O~gb^jVbiy7RGan9TXuHpK?U_leT@@=3h1w@}=JIv+6c`aKB7 zJqXS=+kr5=r@hDC)7u6K&RwTa5fLIyWYY8=$d1U%%&h`ZV*xSR!e?kFjl6PCyQd>s z)Q5Iq<4P*QAHWZj9^b2B?yf+7;G^vNf-Y}|0a7y)2QneI8NPyxcmR17<=Va0cFAAw zyUQCgo6-+k?TEZa5r~^cj~tcWgnPvurMgDDltc2qF$5he+wSC9N*5`lm~0Hru38t-2D~q@*yxN+0$zOMYy2KF_ z4-6$3&iI59TJbsMQN&pke@higGi@%C_>xRa;s$?Az)$mw$F|;TtprgbHh}_`j`VVA zgGFu+!m>&^%N*IoOblt6j0TM?%*!d_5?`h?m*Y@wA4*-CNji}rucx%-zB6xnBB;I; zdcx&ZC!Oj*-=t|WJ9Of+&+$_LDXl%@uxV)Q24uIv#1{JZ$t>jz$^1!~{byA6DQw7U zn~`MRb)dSDaoJ)re_`e@Y222I{4G^OlCA=^$+{MlsvJ9@tnK9S{g9sq2ic~~6Dj0M zUg=o2nUF-izcQa zFqHSW8_FDD#D4;iV#G^yfl_tv0{d+dNdN6qS#iIe5q+S#2}$@ZLg}Lz8Je;h!hJ&; z9c=kBjU6yNg8B~FnHp_I<>&?3y0RjQ0C3{ZFuSrjre;tutZhbGO7 z=g^5L)h~#o>hSMUO)~ANUv}ny1X%Vs3K$;55M+AD@XQ+%vtu}{CQZ^b6W<2h2Y>%} z^!)p^hjhf7`@)MldPwgO1+sH67BK7#{DZ9Ewg&#u;J^_--6Yzt!kT2LZ6Isw&*Cr{tE8hj(*Y~Vc6qRmeBvFuyZgm?$%np-|P@RIB z`wY-^xSJ?}6S+s!8nWF$TkAIX)Od6oZ;+&e+Tb(jlT<-&_3;9!2##-acIgI!8!HYr zLIGEEkjD|MmM?ynXpDC%73sOAiR8I`QUoeXUGgT;oD$vS{?YmzCtgJsqiHP3l!^Eb z1{HJ6Piiv#E!o7m8V9{)9=!M=4NA5BgiOLj?F~Q5QgU_u@y7OwikhXCIE}fbGlnWF z;W}Jl9(oI71fz}#BJmDpcqfv{3dm4Af5{&?R#6-_r%SVy_BX^miTd^+NOE+G-$7Fe zloTCHrWQ&wLY74e1&l#1${K<$k!KXx$s~qy?nBIN)=RT5(wyFy$+Kl;V~pHK8qAxS z>+ck>Dkd0Ul8+1 z2fAr1-?$zK5p4L#>ISj$u_^Z>h+=|bwI)-lNgiwWY|brQppV3)4UFoAi7vJQUmC4g z2$yE+rNOghu_NosbthX0WqK1oWetEna+aMa+4nZZ!r(-wY@?ywomBdFEIYW z&wtW?q4x9txc0#}=ih(zUpDh=Es~k>Lzas;Ob*GP!Htz>IJ<4 zyPB@Su3CdW42u+{=Ap7iX&u-Hyvh=!V8owLA^p8Sqnh;Ku2Jn-Dt7Y&cAm;$fp|qY zRLGaizmFPhRgA2)lczr1!Ttt>lnQsB<;cfI`tDbe6 zf@`nFkv{^ye#k!|3KLr;0YJgCBnnr6OoP9s%t!+!gKB4BHt-9vkRk&A4j!Ivvj%z^ zXf^zK`NY z^FnFk&kQ`sBngp5A(|*d)5uy;e`8U-tC6&Zv6h)$3Ko?kOD13ZCKaULN@pz)^g1+y zH-!IopjO)OO7JU4bOp&X%>G7(KVS^tD*Y*S(#;~{lgLDbysc7Cv+@xB(?Ne>LEaaC zqoA6o5~M>_dgCWZBb0iSiUw6`V6_&#sl><>I^_Ei`L06#FQByAaf*M5yh=lq{FK=` z6+3TH^V?LA<(5iQY9vYRJ@>lDnqCFF4C3%bT&b4bKtZ=h*?*b;mrentV5!xR{j;W1 hb?OUm)ZE%sZN7FKpOa2iWB7ci