From 0e980981379382c712ae6741bcc0ff42f6755d61 Mon Sep 17 00:00:00 2001 From: huizizhang949 Date: Thu, 17 Jul 2025 16:16:03 +0100 Subject: [PATCH 01/18] add functions --- doc/modules/qualitymetrics/example_cutoff.png | Bin 0 -> 27928 bytes doc/modules/qualitymetrics/noise_cutoff.rst | 108 ++++++++++- .../qualitymetrics/misc_metrics.py | 167 +++++++++++++++++- .../qualitymetrics/quality_metric_list.py | 5 + .../tests/test_metrics_functions.py | 16 ++ 5 files changed, 285 insertions(+), 11 deletions(-) create mode 100644 doc/modules/qualitymetrics/example_cutoff.png diff --git a/doc/modules/qualitymetrics/example_cutoff.png b/doc/modules/qualitymetrics/example_cutoff.png new file mode 100644 index 0000000000000000000000000000000000000000..91743d80c2b77937a903f8574f101007f8ac0d00 GIT binary patch literal 27928 zcmce;cRbha`#%1%XGUf=DcK_`TUHranJKa&LiQeoY(k=CrR+o@dq*WBr0f~lo3eiA znoc^^}@2C-V3#aYrX`GP`w`X&Yz1mM(FB-jLT@;^| zx^PsUkb)9@ao4FXV8YN>_TqIcc}4UUVfO#)H#!Pj8c*xqIZcdCP2ovON_u*EjfHci z$ocwK^{4+V_urFfXlS_od;O=*9DKv35^rKGCx^kt#(qHbDxzj@)7VCo;cB_FQsG`t z&SiR08;aS!!dUzx4AAHDEtn*y4TK z`5^pk(^?s992}?Fo*Y$H14F~GL|y~)*8A)>gYSfyn3y_C?Ty=>oF(^I8s_gcjn5QG zvikJ!=Z4Qj^x4i5I|FJOnqyNNy}B1K62b>NYK^?S^n`GeK?)8}l(8R>VWH^r0m>=x# znfJZAmf^Rv*)T9*ocZ#le#OvqijT(ITisj6T&kghT`=EX9^;UeHFdn^?_WbNzLGNHJ*IBdBbTxo~Fb7Mb#oc!|T%Vg=NBta5mV`E{D#X;3FV?)DJ zBVTddq=3ge*Mb9&7IS>VNhs71Aj1IpXqWo zzJ2@lBhwnCNZLeI8@sRhl(54z(F|fo$;iku6hg_y7s^Ixy}Z1<_kM58arKuu3i|J^ zTdH`>k_j-vQFm~7wcQc0fBZ_B1ONQ|eDISe zPY8*Kxb>V^#hn$`=lV3?Ra6lB`1p)XOdN-sJYHN}{LxSLh?sd3VOO<}M|J7$t}lE* z^V^%jx)MKCk!T@0KPCgI8B6F->1ysBQG; ziPY^bHZCe_hK5GN(9kuho%uK8i$j%q{SO@i#XT0E z(hFS=goU7RaBzUwf7`Rw+}CG-OF-~x%WY&XER3+|dOh~o)YJ(SxpjAN)b?Qi@1b*^|1g=O}u%zFE$2mA~F)MppBP8rm2Yx4~ zoP0;OwT{J*Cr)f@RkjeM;OHC)2?+@c4?kb#Fw_b8O=M?#d%JFAgmuJs?NUih*6u{C ztdP^VqVMmYS>MVWSyKY`y{gvQC6mwCU&*QA=H~XC&B+qF`GwPWqt}3tn7Cp0*IJSN zz!{7j9N40p-6t`1bA`1r6D_e5+ba{vQeMUNE0B%MCR^hSOK%-Jsqr)nPN@Z~Z&*Zx zJmk#cH*c^oNrI;DUz;^vx3&(AW;{DSS|5adc(yk$w5m$V%Lgv`)TvY1uH3(dr>3UH zo;cQ=hj3J@Hu(Pi`wX1zoOhCLv%zXf0_p8H5NavjjB`0a?>FFj^Eeq3Q#!>uMk&uMW`{H^@} z?U^%Yyca5#TFYIgG0D?gk4Rfq?o23Hos*%DnNV0@X1&F;B)?VZwRRGMi9lAimU!pa zuU}&=tN~3617(mC^j}OHmOH=jt<26QgoC97$8cr-O^11N1kG(Hr?|e63O5qSBAzqp zp%@ILDfwHs_&7K?P!+2C@PUX~!WDxNbsT1r_SwwW;lnJ~9_*15kh3(sx>hN4qm912 z+GVBuPT$H6NTDc|${b z%;@^h?$-8pWr&ZB--B=on9<0bx9Jq*6==>vB<3P9?!GPM@dPcNv8hC1gqDNV%UU) zC`wC94d2~nuQ}NFk%F^KBU1ylWL5UE+{#?{vzdW6)|yw}-Nu6R(s1txqYxbP@)0Uc zsr0gEcW?+e1O#4v@eB$ID)L&p5oDU$4%M(T-{_6b!d3ul+scmmy@~C*m2?r59tR5r zV#xXA)GrcZVR|16Ph_R~H)JS9F2iY^fC}b4=YO#00a;R^X6ffoZgK&MXA*8aF3JhB z5V))iqBen)$CcAi(hXRyCm8$rGyK~(!>aA^NKg19B!$qN%UR}tNM;7g93KVfng;BN zIlH*<`k6R7I{w<+nCrvBEG{luZ5GVK50D+O~A#mB1Nm#iJNvm$9UWkd+c6z%GJCk8>aBy^VG!WjGc;{M; z@3AWdv=9LKUp5O(sv`hLWw_0KHLUg)dwj+orTJYr`-Df1BzjoSJt*;uR=%8TdTH>T z6Ask%kY~@XJ(_bIt|Cg?PRetv`GpT9^n#7eDLAIEZ>Ou)zo;dP**zSrMKOT68W$g~ zyX8Hr2mKGTvbZpf{&@L( zWY*T!P0h__AL36KzP^5EcdqjR_vPi^bA^pxztUnlAzF9>xbiX}kn#Qg4lgZ+76Lye zR!A1-+l|aPUc?auD2J=r0k8U8{CqU^sg8sE}nVMfvC9;hHv6H`fnF)=R_H ztsg%&KsFhho8v8GT%*g&g^*Ezyy3JwS})QkDIxK(Y-e*J60Py;*A!fuDNQ>|wXxsd zST{B{;y{hVXjNaEuP~)uJZLV2eA#JRat?50PD`q^4=w+dx2METJ$kI42lbx+W3)}6?ef0t~${5p@x_?`ole5av7zItefrl_=3IhIAb^W`PRTO&2h zSsKsc36f{$%^e3y`2)ujUaJ-t71f7P9)In#<-liJbCqTO6$Ec`aB<%$+mA$_nc!ou$-JjKC1Z9puuIVWq{>kVkHsjD5kERuwXc;kNjV-xaD@!40lUIb z-(GlNPd(%}9$i@AKY%lFSZ>mtQKF)cfO5mcY8a^7-AFDQ1(49}YrfG;cNY1Lw&Xid zx<;WgZ$mk|y}Ruc91_AVB&6lv2IYQbDd2z+qM|OZuqNW^Q>H5D6IO4LO5C|~=k=R6 zjR~4E;ZWN+xVQ-Bz7{k>&dFsid}TXa<0si)YL5_1BcOf3I|h~w!9?HSZ|{*YM;N`i zF_s?6a<`y>A5zA>hLBh0cMbT#v6Z*fz zL>&g{KR*_045#8+kEg`&nN-maReH26k2MNA4ogi=PL7To?4O0?vb_39aBicdsOWsH z|K9D%R?6~)kBPji!i5!B^O1rz7_Hq{gqJ33DtI56)cBS+c$#`m(IAL~O42YOCusNh zj7v$y#l=?(-6bIMF&JnRh~!pev14U^pYs3xlPmv92{9fX-Z`i7;K5qRs5Zk@;+WBA z9>ev5wYzuWux2<-ykAwh7fQ~mY-nhh6YM$nmH%p$C*^7D_8vQOBBFZOt^3KzhGV@~ zqCdB{lkw{3%iXwfvch92271V`g0f+GEiDSD0GB5gmIlhmp|xtx)_#uC8X-FhyTNK7 zgQcZ6kSu%q_tA%%`v1a{TM=~GS{zEKm4sZ85X(hIOe_yc%*1y^DKjhUIUL<>z!0~W zN6+6HD800l+LP;>ii3?!K}ToZIfj)<8bf{HyK<|U|;gXtD~cX?qF{zVCo&TxPs8e9*47g zQoy9@R4~Oj1S4uV8UR?s)^xo#Q|Kx*m0jsIFg>>HpuZS>c9=pyKmZ~!tEVD?OEUsK zjhKLdU~F=dQ#qPJ??<9G)&p;Q)VjI|_+q?WB1za$S26SbP1gAplL}XrCGe z@FL;cqjGWg?`!3z78Vu`Pbb0_IY@1_tiy+eLpK@*Z%D+*$oPvmC)ID$=#jAfx7X*t z_4h|SctBmfS!B_IR%3Cn{PLLH7eCYe{T;_@_aPrh42EM9YG2^qf3zK=q@ZAXd!^vE zn;VT_&2~m>TU*bZ@vzhqz@Ksx*&&~-%Iy9S61pZhVcd_;Itrklb$6#FmtMdBt;*|7 zy_3G;mj8p2nf>>!A|fJ9YJ%0AQqs~;44)zRT%MYcK4KjZ5CEX!=7T=@_ss}*DxH(j zje>yA1KpIoIZOG-7cKAhqhbG z%M^uGs}aPM+HNnhvPj_mi*9^oC{bX?x$P^&A44vN%9YU8${oBgA~-qf7>YmWI#NxCi*76g_m^Y83;C~_8E?kj|VG1rwrI#2{qY)EE4?uxkZj?L2DIYU(M~` zzas_*ucEhiI(++Hczbu@t(Egf8<@yRkMo#9Le+nDH5`!$GG4KyZy@_mPERLQw8j_F z(9r04D}DQ@QOwDa%M(Kmuf6l22|G9}O2ebUTq)T*ENqq~h8*`bhMj}zo$oK_oH@F| z!pSyaMMd4Sw|&_H@d*E-t+KuQqPPH1b8zJ1)3WL2w6?bD&?)EULJw25+&}_2_ijmv z7>n=nc}yK#HDY1jRaPp*b3!oXzF&E0-4erm2<99e$pN&rrARYAdh{q!W^rpIpoO=v z3LCaa)TR>y$k}syEOaW8C~^Js*%5?cp&#$G5xwo<(GrM-LviewQSMSh2uarqH9Ba0 z8-UY6c@VbRY3d`*^Qz(5+?h3p{%9z4kG;}RD?R#>|) z;k_}7clPXAI34xSS3*A|=+<>EqvC!1`$$wufB{=xTMOb;(W@!mT&g{Yfg{k++>8el zz)aV3f*JtrB0I$Fuc4QeTJKcIHYj9y`t&JcwS0d~C3)|aIt)b+QJqYKH1@i77^zSt`9TDAUu)tgYG^RJd4Rk%lbe&T)f7+}d# z`vG$R1MhscUYU&K?xxBeXukMf)_OPi-o1wv6nc6ES4vbd6mWSpP~?-(x%`NtkPidg zLjq(V>UW)>^sNjy)*d@PsRQqI5)+1+a41JBW1+0X0OKg>=?#BcjKIov0VF80{BYE1 zd)bWg*fH+AK|!CMs=B~tn?nLLhyQb&Iio*6Ty0FsB4vo7NcI1nV>0ryJ1eqm#NQb8 z*&HWND*FP;Z}{{{&C}Bp@iz<`8wN^BsapX0KURH< zLYe?Rl6vm~4HW`s2cRw;P$%fV;YZOWRr@bhZ>8)qi+e7oZBHg>1}7yYEpKellsS$( z-0|Dn=xbA_JTtj51vr;5$maTSX`6g2@V@rvt0zn;|A zpi2f+wig0&hH*_ziet00Jcz>ql*A5w1sAkWCplG)nnG#G+3$XNiHadw8Awb`-?h(& zJ^wfqLv)h7_kbzolyIG41FXS$<*kkO?AM2r2 zHYNQ!3jj%v+7FZ#bd~a6*3x1_1Pibb8*up}0+=v3$;7M&a10{>g6S+Eehb3Y%EYZK z6At%E9S(UIKxBOeJ4?$<4BMGo-?VO;Dp*-r>E+Z4%Q#w;d#$}h#SJd2*Uq%qGW3Z8X}XI?GFb}U zT_jwcPR9$M5kQbvUi4kmePZ~x>p^jWHCAj3g@Y@O$3YI|wq!Xt#!|j>)aK zda?R0V+D>6@o{Wp;{lzH zjZOt}EDUVajpe!HA_Q)3Rj-S4_iN=KZAc^VPQc+E10X|ZzfZDBjExPmv$fO$ zV#2zN!LwqrdiWG3 zCU8b~fEfI=U}9oIA%DFmTN_AI!XF)}v2b|C+a8@Zd}~7)!zzo|LxtyEqX6O%ts4ZG zANpw=%NwO>WteZ?yonkFOd6bOOhjxoHNJZ14LEN%8)V>yFllhoycdAjvFUyp zT~(6{yniQ9r%-K^u4Gkq>Hu#c^Lys_?*7h_=gvaer$LSL=dqyXfogz<$}VL8hGb!j zTpB7Z-qMp zTE_v;gs3S{Shykkg<6;T@tE^TKxP+?`tz?(^AwEuTJ-01cez?+&D~!gJ;5(Vsta z`d;k_%&?>;B@tktsJ|#p@3OhQGOh#}p$?k&!Rd&Yn0wLDN3N9E>Y1r4{s^8>+L)Z0 zs_*SB(XA8+MQt+TT;Pt-TB2*6`k0955r{KDF+4<7n6Ujk`~ss8G1g7 z3Mt+k4hpdE^+4pJRV#Od#{B#WoXAg036Qfna&mH#rG3r-pu7j993yp5%%abC2msX` z{q+F(XIaDTbRvmoM|>}1D73Y;r#8bOhhr4g)uT*&%*`_(yW012mL-;{?LmIV4m5;z z%wz{xM9^9M_q@vs@`s+t6C$wR1@Eb|P*YQbE%)LRFee}Vq9xQG^1*vKJ3HsN1rdeZV!ZxYBCD%|EC0b?|o8Dg9oB*U(^`0Dh1Sp#6G3&^5f6)=+|<&Q%Es2_P31eJi$# zgo765ITk{v$~^2PYr(uZ34~96K>@F{^-JJA00~Dz)giXCv#WL);fa{q-*Sf*96;X# z2_~^qr!Jm5cP^~bV`*0U>-`s#9&MdfCUJn;A|T;?3cZxL^{jHF#I7$mG&B_W+9$(f z&Q>||2f-CFK619{K$xci>;nz-zDnwaE#A7l(8(yqvt%98TTHNwceyXW)U_MO)V|Q&z_~L-T%!CI4EX+Z|C@#Go^#oWQk5rPW1qgRZ9*auss3aP%<+| z_w?viZ_NVX%wH1PO3~HTg`~kMwJThloVdVKD7d=5t-6pVk;<*LGHuanZQVL??b@}K zp39nM@UMsuhcFI>X4?!DBGB!W&CJZoYoCHvNUSNfg6q9GAHKZ2tUAFi!Ug#~JLlNx z8>62DYnx7YNCh42Z`GnK0yr=TdLFFpX)T_i^z5HMf1Y!lp@5eXVL&Xx=;-LYcAip& z(plR>s{Z~f99ogB_Wkv4O_VW#xZcqTCE9~k0RVj=F_l+Bf&rF2I5sxc#0Ml`GNyBZ zkWwu`KE?!YZ+ixn{+yo9@~AyI(Eur78IC^UqI17&2kdOxBauU%a30$6hU%O5c2l+b ziZD4Ls4W$}`K4*H!v$E0Z{NPPPV>)PjK|f~)C9gRuYbS+)~vp=#|^sAS7$>O;Ji`L z(w<~(^nxD;$T}mv4>{0kQCYA95)x7qaPjfYAy%hO7D6WYSjOrzFN}+a$1W;L4VT&% z`ThNS<(lz1#6Phs+^g&8Api7`zY)^Yk(j*9V&bH>Hud#|0Xk?wJm(6^I&<}SAx7D_ zxrr>^M*)s-FG+#TF}Bj6%3+*Q^C!q%djzPA0{9uwd|5z1U;mAmTo8X~9P)rpJ?Hn^ zHD?pp3Y}1`eB&}QdBVdiS+$G}f42+^jrdNGdwJFjM**J`!I1@x@C3*4c&(=7o z5ZoKbV~P^R1g=(l(;&$PE_!NB7x?BC$iYY)RGyGG>M;2FF71hVxd_Ops|&f48#`7v zZoofM12Ys{R#xV_AGatp2sNgl#Y<31Bne4Lo8VA|Ku019LhInP&^kJiTR}we>jCeL()q+mn)8 z@c`&ouDgR247nvqt!Eu?t2zWuoWMKlIGn%!X6{?jaTJ}InVB;sbEO8WD}b~bA)X_k zOBx4i6yEjhYfxu?dH&ep!runHTgax94H$!QAOe9Uc_l>4;i7BMLSjKkMy3pR?*D6w z4)Jc&9jOFbPlF&xNiaZcd&^!-%;bHn?lk^u?;``<=HgIfem)<-F2;nDm*i7r{5m0p z@wvXQG?>&d%s&J^D*Z=5nNdRL>jf{*yldHe1Dg9zA84&sP z&b0Fv)zuYsb#+57kn~K#X5B;S6I6rUosk1&S6A0Zq2rkQ0VOgnKS){p*0sy4^CM6F z^zQ{K@E$z{k*^yM4XHRKd+XXCxq_SOXYf@U(;g z%C^|uUUA+*pbQREkjA!*@3s(>k`N$qt^CV^C*ysc_Ejz~QEP7NCuOo2#zwT)gIQT!#!W9M$O0h9}ZH<+uQHZF4r0WohlC2f}J zuHLh%R}&z!0TrjS0NF?3>QyF4rD+CiKO;THI zs@wqQ5yi^N3dh;A$j1|oRv6KRuU`vM($JjH&o^X8geUAq#@Qt|Jxw0bQcc-CT?_?k zdQp21_nHCv)C&PcdPYW6L6D=DL-@X+sUjjYv&& zl<@PwUQLov`UBhGTnjC<5h=Zp9H2gZCMi!XbTj#eMe;zb8;aWYT)bN2Ya}QDkTOU^ znQv*Y3ff!-Kx?IIuxAHxwpLrmdZ(K?U(4@%f6D?4qczMHK)6D2nY55Gb ze&FTRY>$(y4ErpFXJ?;UfnpX3dDRF?R2ZF!Ctxq&MHQiQ3CYXL+Xf{~&5u>*^5y8& zitqI-RM0=28 z8@jr*-`e)7A}(H0Svdl9*W5N0!rV`tRW+U~=4STx{N@%Gyn`l2Bjr9@i{XHhpfhiR zx~VqxJDzT;;Nt<2ZWP`{ez)m%W?{RpeCQU;?Cf}oqzq`Zhn^_tR~x>6 zH&TBrs`3<^96f>YPoD6e)J)++_XJE1MiCG>@M81xw@?@pK_sB}LgWcTyujo>r;F{! z3!^|DzA%kv5=6kB0m;m5;{Ef5>Zt8$R5R)ZN0%s$N*e`$Hq91fV+v%vs3(1PU)T zXeZgApQaA4L3$5UNjMn}?gA2M-jN9gkq{;TYDEMsr@;^fXvyfOUH;6T6cBbZy?IH< zVuQRz!0*$79GClL+Ajc$Blyt01B3!K`T%GV5%fmPQ;zGs`?=62m4N>v5$Fjj5S1uF zxa4sYmK&0lIbj+}(M)&u>RX#Qn=hH+=m(&EI12okNzU}QH`cL;xmBw6-U!yuqZz7j zqhdUB>y(=l_t4sdvcZF0#UHO;#n4so>KSZ$T-A3d2aHdTj7B}#ap3-=L7xTHoH_t3 z0FOG{HfHq@JB;!KiV?W_Cdhtj8+kTW+~NXOmX?GVa{%axHSMd!6zst$u9cwbhg4Tr zdyd{AXafC>2U3-m?JeMHOFeg2T97SmvO{lr{~M?yK%JO3g;Aoth4LH@aDh?mEW_hI zx25U5)9;1~^Ov4FJSurQvNjK{00 ztIe$RHujFW6=RS=6PzMFWInlq?!Yl4xkGM)lWF>u(I;aVL2S)O5(}CcXZYF{CJ!roO{SxR;{dYu&|od z+)MkrDB|LTt#bpt5-)1sFbage%KY1FGt`}-_L>7a4I$6v>*be#?+2GA;=aF_#or00 zxY(*LCw~Jwn8G&H;#WjO;JbI@FS0K#*Z?(}d|I!rN;UJ;8%Y97AS%zf%~Fp~OsEy8 zfrO7i*uTWb({mIA)hUP5(3X{#m#^-vue%;sQbWdGMLjcpA2}emgdiJ1S~~)eermJw zAefF_F@&ogn#x89QLw;E(AsxGW_B*i>#~{rNiAN$K^bQSE<^W(r}r}gOr6(Z0aZT( zY`caPuH4X(AjU;Z{v-!S2-2+)-hkwT!5Ej`st07V(yb|*4jsa(=mOwcf0$iQ(#qWU5e(yEwjua(u@Z3*G2m-id(UBsJI&WY8>|Bj*vfH`t92_r-Fzgi_fj26xEgtW+0UJ$hpvA0H+l`VW-qFO2~Bq4$vEMyJ#DtMe7Vn@ zBb`^iyhvJxU-yk482`e_HWwGKdu`UNrV#_@f%_w9>7dxl0br$IVTplkXaOw*fD~M< zjhp>gYfo~YiC@NR%TtUPXQ2=&=gaI<4ynKX)7CSEPVDMyaV=zVQ=g-KWDs`!H zk%?kRbiuK9&pl)S=sWmNHQ0{ghdPZ=;Zg1gla&f^SI0k;O24J2cMGyvdzonuyBBo_iRJ5c- zdc&AE5b(6L_W&hYVMCPLj;YKX9Hqc;S~g4Pa;5U*o$0G>doBX7roE6?rY12H?`Z_Z#ceMrC@Q~y znuGmAk3$JIErIu_(mCfX!k9tnyByk??&lsBzNP@L+W3ptTK|zW{`m_@>d7m@QM7rcT8wZtHyDrum zPh5s4AuY-$l-%#zaaT3opvmh*gwUr+#3wvp1_by!y|yY7qG61UG^v=(IYi0vH-|{6b9AueIqZ86-kF zI~P|0e9NU5)_bW%Z9cN#_szqHl~^ggf6TDvE;byP@fT@w|FKXgyF$(O-b&&#x(gRP z2C}QoNR2c`?s->TclYJ!Wi&7mmW|)$Kz>mL4kqaH1jy9}`gytZnKQIt{zJ}Iz^HeD z@uGRo`KQ1Pu~1;3jpswR4f;$^p((4+=KT3e_XP~*{Q2|rXC0V){ruvRlf%vqd*Xnx zOio@NTa@9_qk&JKKPx$eTS#Q$n8{uT zX~>}0ZkP%R3Zm|Pyg72p!Et~07f37#Ntv0PmqRoTUrZzPNt$hoN`8h}UhfO^ThnpoH)B`pUX}_Aa`~JB;AT4iny+SDyu;tG8{r*CK{iR zfU+j!#4#YOftFE( zvh9WQiYFWX9IZ#!u!@1JaqhcWa)n94^A>ZRm` zO0~1gozBkj7f)TDhqXHYeR+Ixt&_w`>@u~##L=&R#oa>x{}6Y$XvMNca9)}wI$kKW z&oJjy)YTir`M;dYp=FnnVglP00J!7OWC+{#oI-uQw;y;9{ytra94DbJ-=^gg5syy9 zH!VB30%tn27dt#0W!K|Sg((;qB4Dwer=?*?=md~7y8;|)8WdpXTHpt-qt4k6nhUXp zf;0b~0TM0o9HYxnlGE3h>wXDqnfVU#?Cj9U+nZEyIpEs~yujk6=+I{)3K?HlLhAW{ zWC){7D0i4Z=p9Xbz%320=mM%y{gDVB$QU#ZT!3eQto0_LqxoB!#BU8lsNBI zfxN$V4LZ$LQlq^Fkx5ljCC}M)LMX6ZezcuOWKp^TLAf;401zxMxVnl$TYm?7QIshF zAWb>&!m5F?y$VR`d6U{Hx<8ru1QC01|M^q=($dWEuMqh`K=SvG^>}{|`98De!iuT$ zcvc(qwZAI1P;4h}?uHfC#vI*9dmFAQ(~gEE-;We71qB6k?0mO|JP-;l-$Ma9Hy6n) z>5hhztTsP<{@lx}g^E3VbibW}r1fJzvfTjX=mf9CY*kwc}7oRD3}<%c@-9Rx*1MPUEPTNj3M zE3cz-6q-fqu=`6CA+xVMgcYNda5|Qw4@Y zGcbh0xoH5p3(ZQvaElq3y?8D_@{suJ^ALxNp`f5h>#G4DNpNs*9hhPPZ!|-0a@tzF z2J>9$0SEiv%3ZjCX4DHy4h#F(0G6LHK6xvbmo|SLu5>tD*%jxG2>w+mfRAgP6xC`clwghM4HZ;>V9pCxNxCAh+FFL zEQN(pH%dc2{iXElTXHcaX8Gx~%j=S2J#h{a(ERwWUGep)$F6JF81f8j&Hfu*qXLYV zftkDf6$g{D;)!26di~1rU%%?co14o_{St9bX!AGzog{|X4S`5jfVOdTauSz|i)$Uk zXW$3#_mQvzNWcek5!@m!F77ByZ}2FqAEwQFB9+*gU_nHa7s#&+ZL@h@AQrv2bLLe` zaI7-%8oX+Po|}S|l^ke);O~L`R+wr9^9@Q5Hw*%!=!LN%xhnxUMT9(`Q86x%wk==? zV+%(9oY&wDhjfn*3~?P`-`h*W((uXI%o2y|tJeUPI3C=%pe&JDT3MyT;1FU&fsw^0 zXW`rGvJsu>D`W+*^wK=@@QN)7-vBEDL4_Jk5cU4J0tk?iF_M&=oVKuMDI|cPt*=+E5_U$su0sim5BzFk4DZA1$_1D(<+mA&gy;nFdh37s4^8f zxRYe?&D(;YF3T>5dicL5GKMI$iOa-~>#$9n;x zp4X-L6&?nX#2g}&e{$kSSmOT!_C(S*=#5H%*Sm_X)Ijj22AfOrIj4;N^uLzuMB6J& z`TvSaF~DRCz(9{Ms42+j4Z3OFSYzne%uKfPq;3i%N8imYYfK;)qco8*uSsv2^6F|} zZ(gKh$Kz(+&7-{{ZI5`ao1D;j@STPY=#%W*)toAXBqR!khRGsW+R-#&*!lZ;=Za!o z7zvDT3l+re`d(ig93SY;a!KRNwfg&{ys07ki(2RvxKTzB+!}AxlP`z@hG&!oG8cU5 zkdplU4S7BMIoB%6?VOEX+2!KeE1swL;(ks1^MhM(q`w<|AAA0g^+Q?Cz4?VZ$;kl@ zRN4i9S7Zcs-#pCk`dst`m6@+!fzbPpGtVcNtx?ou5yUTi13oHa6>NkFpXh#lb8#WZ zTz@-Y^G&rX*D-o>8fUJLZX`v5=JnK$itdxYGsa=5jv!uk|1y9HqmxP)al20BgpZy^TL#s%`*rn|3*O@V zS|xVF0*91Qzv|oaUS`Gd|E1aoQ@}m zj^fq{79A~=95)x5%Tm1~4CY@*OaNx^-a2dFkL%&#!O6)v8h3%^4h1&kss`Xyg+Yqcr^-M>s0<4X$K|( z)!a8QkcmOoObpDaD8kvl4z?m^)@bD*|GnR10LG_i*`TsRpP0Cmm6i2GX!hg3;8<9j z;)k35&->2g^aq{znK~d(&&|QHb|SbA@3{Gd z-c^c<>*R;LT%-sTX=FG=HYtX)4t!-k-rm@!hN;G38nG+}L$3fWxvzPbW-aomkoL^9o96I)yp~|u^ zV@D~q#UZ5#TRv98(2H0S0vMS(AtIvlA}l=YB2jvJ`nOVhMhNB6nVB%~OT)6VNlG%x zt)v@(gZDj%FKCS8siX%PJ3Bk@y>$>edNCS#fWRxKpxH z!fa2te125Ifcn_B3QX6fH`_xVBR@-GS!GRXh5%TimI>bYTWvtC@5aF4Zi~D8Xx)NF z4M3emMcrRDTm6P{;DfWD&5ZzqjRFJ&2lm7Al*@d%d2-ZjKO9IgX#qj|ujAS_?W~Ux z69>Yqw__~JwcSt5l)Vfk!q$)1@F=hlsnIsS-w;fZEbej=DQTd*=ydZvoPzbfsH*zq z$=0*B3;xbM8-Ccv-aZn z@-Q(1a#L@A{}DgZm+r6*%m-e2vR{5lv4r?<&ZJHGq&4sN@uGXNIncV9*^vYXACf9z z0@h^7{xo^_N$W(XOY}Cbo-S#>ba|5FUPTp;O5K|d+A=*6ZZb+e*27W+ch*c_fVufK z{}}L+dbxYh_1JvPe^iumk5B8)sk*Y8BjW>a9ktnr5dlF%4R$^{K|EjC54%pj*$mg( z+qT4s*lBW^HHOkYFDRuUFKDsO6p2$W`H;*zcIYC5Zq6 zMmww6-MmYe`_YAN{@SAW^HbpFb8DQmIbdOX+W=XjnWT9!#^>RtBN*!1z;D#R#D*SR&w&fB>eqbAbT_H5fqP$cT^H>{f*VD5OT9 zp+;Z|?=JY&w!u*H7eh@YS|FLmv)WQIm}|2hcS@S_et!+t(`J`Ac|h2-Pqe;3Vz{=) zx%Pe~-PTq7O&u@IS8jWDqW$P6SIyX$7s{Y%*0qM1gcL&h*YE<2U#x%a?WH2iF*MO; zd3o(AK4$yXng@_=dNHlV-)Vk~4;7xom}0fI;m-D6mgLMJM!_j=O@{OIOaMOhYm9LD zQ>L_l0AS>DPs8TGU?;uZWw@KjKX>E)=KyM3MU@g{LZe2rs2Up+bLQ)1$~cX^-@QIe z(|KHYPJi9He&9^ZQtyY3juv3>(JVFC?-=K9LfhH!=k4?Zt9BMEtx7{-R_*lZZbMB- z0JQfe6h0?bWq7EZC|tZ|!g~4goYy55oGogqM^PU-upWutJc`2T_5qYw*QH@*uByRa>Rs!|8sUNq`t9%<>d5I}G;tD6LZS)B0P%9iuGl-ZguLf&046fhY9tWAv$?vn?_@sP+-iK-N{&_vONNT zEP-3O05l$p*CZv(hS~P=@PWzSza5nJIxspx&P4P6pG`@RB|>$Bc2jvl$Tk|=PI7*J z75n8gLIX)^jv#tF-)vr2ahElt4RZV3)ARFz2BvsF+Vv!lW=PcA_{cqd_H2XJz`y|d zOn_=fLCwhBbnpIuK}l8)$h!B3Q{iq*o@fM18#2Meb3T9=n;)t)aP<2Z;JDdj0DCEB z^;VhJcYyPx?j|Pay$yC(SlEp|4j$j zS=Lj!A>r9{6yjcM=ew z_$ze5w_&UckuWd54#{R!QDh8y{Q$#J?7!9N%M0(aZGq8SZgXr7z3vC4)25(==Yr+; zd199$hdvA2(DS-h@99~=s}C%ut{GrW98JxR*?Tz>;?^Y{1}F8*NUbcIWJQy2H+}w%*Fh7B&Xo>YpqSo)$&lf-ky-=A3x}An zu#!@OgARrr5WGZFkRYAnlhOx*WzCtjaSb;->;Vc0^ymU`tp$Na^rJm+0*wT$zG)r} zo>uY-3J~FUzyp00PW7rdD{Ld2V({nv0hOrvg-6XCh0Ovx5Dq-oz$QNl;!N2*f2E!R z{E1199)QR}JREJbJtM%ZCZpMTpJe$_V+eNii<5;fGFtklM-(#C8OOcc6~%76 z#m)pN8w>|8`t9xQwRLsj*8jg?(mZBa$S(Ja-iOYjX$AOTu%a6}u%>WbRLtaqIxlv; z4%7Bnj0Dv3qyaR9j)rv_K%YD^97B=OrexrYY352pIS_IJI3dF#BNMDYK$TpB+e()J zxDf&@dbJDwQdXA@5J#D1P`!WB+mp8f>w}(qAZp)F4U!eC6HqfsK$oLIC=2N$z1k%| zDX?L}3l0i+|3UiW9MJ?&3+Wg}fQTUl{I8tddw(NH`7taod=H@E%t^UCHK2PMPyfufFZO5dYQEAr%t0=`=e$-qe{U@bj`T50(0~w$ zv=QKhzE!yKf<{OJ#eRKhgcYLhZMK4|o4AXbC4*KIw`U^Ktn(QT&eDSg7-G*nobioV|1A*GZ4_`mhc2{Kqx+~DWZXe1&bkKY!zAc6OeJDu2 zVZ<7{(l8GbJ)!?Zz0V1Ca)ql>tT(jk>d-NhfJM zgf%)7*ptYr2wCOhHCYhgp+W6gzu9`1CY)n8ZjWHUbstBm;qB?+9=!Nd^4$ZqdCd6Q zFJM43R=}F}%sXPNAL=I%63=^&*9x|!k55gsTi(CNM{bf!srUE{9kc~fwl|2RR_3T$ z`FRe^uUpbzI!|_HkyIBdi7FIl7C?@^zO>g90LQptLt{~F>r2*(26W&^fSH%p>Xj_? zh9TXo(6KEv!NVE$vv#Kjh_2Zku%urR%QO2BdGu&py=@hNB^VO6w`vcfVQ5ch4#=P) zn=TGgR+*}onhLrnkcdF{J0O0&Lp~D|Lk@fYkYNeS3!KL0k`-YY7{8CaXOMrRGB-ix;Ed0eK!^SZdcz5Q}8xJ=X)lTI(RN8v`~ z_mX94AP$QbqXNjsIx9MjxQJvsVj!4J7%&~=_Zsg83C^vi#) z6H>r`aANli1}bU0xmqVfld-bkzUHdY+^4+L8iD(Kb5~^0mipLC+GOfrwJ&YjUAzIGTfd=@%{*Aw#G>Hj)!NVc_tQZ59~_X2s7Ls*#7v1$b$ z)wSSYTpQ8HBCJ&4PJl-X@iKAi5o8S`gHY}e-z&^9P*r^^{|8mI1SemG7hl2>kV8lG$ihCOgP7Y)iVM86!YK906mF_S4$iXbx+tg!50rw9GhHuNv4 zFHxvcnhge$-?{#b_mzd~1NY3mygc5k*RD~5^8k(U!ZZ_lloRNxqd$J!gPqQ!7egI5 zbs1JxE1)^@M4D~XNfdoHu`v+L_ss){tg_8qPSdgkCs;NIQM&J>fp;gzb2 zA?3uKUt7vSUWw$_AMRjWUUD4zC9s3l1+q>N`1vme$;fV<*U*lEsFm8edrQ&(i#{9# zU8nUKZ6uFCzX0VEyjs1+g28g|cqDMG!OSPj_XmPJkVkdv$cs=ZqHmMvVo-tNJ9Kym zBOMq6`H4qNted%ZYnU2f{#S&oDK&EPwe;4cR@VUO#y-l1}YNE zH3;vZv1_QR^rvr-VvrdJNJMDKaiRTwX_*Z_ONQ9Lj4SN0r+}A{FLCf?mrW-OQ87un z=MT0_{rJHlEq%=q3RX9)iQ>WNgDT5<@QPtDkcGV$;J4L^;=iq0P|4L9xZ>)Q6LKCN zE@erC>xit~p_=j;A`%eiYW#x-ma`P$VY8y(MNw2NxeIVFa9o-l!p+&o*48$sn&ZSx zBqA~-v|LbBY(}}TUA4HO4j)JOb-kRS%I?Qr-E_;4cQ+K^FtyGDY^;Zc07tg7@jt#T zN1Ee#EKg*8u7GA>rWj72%$ehk9|V}PCgv@IUn?pg8*3?=!p$P0p78_w@Zwod1_^z&)gEe?s$>`pjc3 z2U*?^#N?9FQLOD#*JobVveESaEMNsO?1OoJH1P@Er#4?$r3>F@GZRxjy8l@(KyKMX z{t%5sK8jY?uy|4L!#}K0-B_5!Ba=P5YCm6e{o4)5n;^2;0fx~Qe&|-gB8G-8;M?oS|na~;ZS4;In;y(D@A z$CTP$(Dtz*Defh^6Q=(8jKKoo>;vMp+<1dC(|SEQWGHRYsQgRGEQ7#CjS>K8-Yq8} z;2MFQj>00(65HVWMCE8)00KpD_UVMs}6OX#Fh+V z3g(|RLka)m*@Ab>CIv3@q7qa$Eb1kqmDSbwfy0QbAN=x ziaVFfIaTfk_yrhpl@Z`tZ)N9jd?3!M_Q_?&r+h^3^3wS)To8KwnUvR2cr;+rQU~?V2{5bG@oO@sfzQTU_#C*YgQBAKmsiB+Rf|hNF6r)eu zv2v=eIp1EgkTCV^{r$2Yvh;QY3LXL|z{0|!kbh7K(k^;t5g(s=dD`Kkx1NHyWB+8} z%rhMPiik&9Ub&IP0>xuo@{l;2?v5qh4U5r7cmAK=&O9FLbp7KGrBK-_9)?2bq!2Af zDI9xtLiT-Bq{47)*|+SeMkr0Tcv7USqg9rSMoo{JsEmDD6sDRnD(7=O({$=}=J$I2 z&OgUrPp``JUGDq3ulu?_pZ9GPt;nP~_Q%Y7{+O*#8O{x4Z!Hw45kmDJDm5(F)stED zW>e+18;N}x%;wrlLRz`f{t#@wovR@?P+mMX7*S`xx6 z-Dds3BGWPUm za%Eu!g?^+Zfp75yc1(6Z!tnWLxrih<$cZ+Dzxp|r)K#eUlhx-?PEH)0szcapk7z1< zqN{Pgp~pX@&v$bnqVL_qF+Oi@xd*q;7cvdqATK>oV`73JrDj4BI2elH)=q(`#sg9< zP8c{ca|{)~6+o-HLc4r~znF+~ks4-_zuw0=eVjgdO^0r^vD&a=JyI}3N5}kT z0*o+ao*0Q|Fs?nS!)TCcJ%7H}1416$4d*S1y2UU9k2-yLdx=NJE(UEa(aGrbuy1yD zQc^Jg-SxBJfPC8(#4r@HD96%zJ27Wj^xFW}5zZ>sp=1AW6WuAeEQ?`S(SZ5XDN zC+r4iv@VsdeWF+&Y&P1l!)jl@)0UENlAiMiR{OqG|7;aq^FnA%jdNm#GFOhJ#=^9n z@~O8+U#>%uP)`BETuGXZ4iA?JvFy#tx|Q~II-8f8C|EsCtkl~yyZ`8OH`nE|MZKNb zX|p}U+_7y<4bnjy#A80{0zUS^{_2B~s`bVy?R{(hIgXtTmOumv2UxXKR@77u6Gs`)+^k7bI-fU#k|L))-G(ND*v+Wkd>7ve0nt_!+z|{Z|-%}-C-es z}Kd$=w?esdGvg1oe@P4-o@L~#gnv?bf^vFW?!q_=S;!z$)LjEMF~epqn-t?*7$(!O`j zWr3EK7Wf}61OkqlIom6PXOu?<81Kt#Y8Ve7lk@2=%z))&5YAei{>ulb0t*lcOG-jm zxFVhWEnycrxV=G|2z@_2MmAm|+Cp6smz5=m$a^|u7oMV8rm&vg=&eEiY?O6T&)}!|J(`^CjEsebg7KZH^->z|*Y3JzjcJ$|fbJsL154TJB6_ll1)9W4? z8j6R06eKM1RnlJ``;KW9)Rkvk@k?0kma>le(XU&MM!MDxs#q)9|=DqiLQksWJsw2|@9jITtb zTSgV@Y+gFL2Vn6jsNevYrtq%6^^p z7bY|#sX%njw}6vQyfRQs;pwMhCRi&KQu(b8s53&62UD{k9$h5di@wR;cVIjwq&Q@& zKz~cY2Y($-kg^c_T5}8@xoy|6nd}iLHMg(v=KF&~lQ=0qnwp}r6v2jllRZU(mBhf>IfX>u{(?_N@2S5PDUXi8pGZ`1+WnGx z6k-;4Fu&3C!GQ^R#|OOkrLj=Ir>FTj zxvHm!7uU9|UA3W&_XStQ5HwkUO|wFnqAQ>}0EAgnty4) zU<-Wb&6ZfvFuz2U)N^v%mMWu{%Z})RtX35;UlI_ZsAjh+@H=@JRSwALRcvGxQZoZR zP{MYy1zHE?<>d+ANuJ<_SUfl&yt#s9vWp z?IV0Wu{H-#=YYukJe6=!#4IT*qk|jvu0}y0wLvOi4Udjq@9lMLO3TVGT&?Nf^0{U0)d+}@AW$fT1f5#60R5e z+SN6K=G82ibK|dE(L(=lZTPEQjaxCkVZX6aRr^jM;VQ0;L}NGF+Bjn|+cy#LPqM^@ z`$HzT1|ZRRun00*4L8`x_uupBN~~S}OQ5RET#%8ms>rb=SJ=m}Wzb|Z1ie0g{q8vO zd;bP3@j2(@~zDvQRNELaPmw7QNOE(_&1?|4Zr1lFHv!zS-bdWTBGU7c`$3ozCUH zn!o)thxiF(IEE7X*~5zoZFY};4a8IRO}2ZsU&3_+!iwDKOPL1`T4}5=t}tm_tDsQU z)xAz#-G2PGDFp$-Yrfa~JJ0|N;XalGE&#!*l+BcZZDF?w+`{-N^4!7|!->F4{#~Z} zSg-moKN#slXVeO9`S8piQ8RUm7;;CIOhKCh?1_JmMv)~R(9NWU;zf=N;cMsR=jUcP zgCP#&92OB)E84mc9;IjS@ZOp2YY153^ZoIIFE7sspv7q`KMz88)F5pO3!pt}n~d#f zK;FnjEHk`xEsBRsP@uuXUsIKqk&!63jj%l+GjPMIfCwlJk~C!Z+CXxEESmr#r3j&R zmy4m?OVbKxQGl(%C9!c&%6+$JDV>zSp}#ezJp$~80nW1+3Sg*whkc6{+`Wj0P@V{h zh{*PXqR+RjvOhoUhD7@ ztQjR{Se*5si)acNC-N>0@_AZ6&EoO!0R*>wTsqr7 z`(6x{`&J^RfWRd2HWGzv@%H(dg5Hs+sRlpBN}mLiwlE0T35i_f{drFwXtFEo_rifg06rL-Eyu@C!XdzZ2lX^!dD9Q=y$>F6 zX;-+lfD(}a2;g*-hldQVA{ihqoWcRpHz@DmNr(`$F=o(h5^eg`T4f+;0s0;bupKb| zjwf;rQ1GGSOSH>S#4b07Pii@ef$=E%EQ5jMwiE53_xxCW%An7}*!Vgnnb2hofCS%y zSeylHxcDo3NWIC+vn5mr6$EWZM+t^IEG;y!!K2QClgu5s?C^L)+#0Q?^GxQZmb-XiN*J0S=@ej?R>J#fOytpx_6J2}oppo5# zP0b*J3>@(`we_<&W`o4bbk|UzsG^-7Z}HOw%qarrfpo!((iYq{atu1zLU;fEwS{4$ z;#4Sx6mqJ{X)_`83NGy$ z=h8Ce6!X8zPhly6|0syloxl$$MZZkX+dy5@(xv@B?Y_&XhlHQhfsamzkImr!;fJnS aF4e7uJT=Gb`EW=XW1GIY-Z>rjxPJjt)lIVi literal 0 HcmV?d00001 diff --git a/doc/modules/qualitymetrics/noise_cutoff.rst b/doc/modules/qualitymetrics/noise_cutoff.rst index 429c0b94b3..a96bd4e6c2 100644 --- a/doc/modules/qualitymetrics/noise_cutoff.rst +++ b/doc/modules/qualitymetrics/noise_cutoff.rst @@ -1,30 +1,124 @@ -Noise cutoff (not currently implemented) +Noise cutoff ======================================== Calculation ----------- -Metric describing whether an amplitude distribution is cut off, similar to _amp_cutoff :ref:`amplitude cutoff ` but without a Gaussian assumption. -A histogram of amplitudes is created and quantifies the distance between the low tail, mean number of spikes and high tail in terms of standard deviations. +Metric describing whether an amplitude distribution is cut off as it approaches zero, similar to _amp_cutoff :ref:`amplitude cutoff ` but without a Gaussian assumption. -A SpikeInterface implementation is not yet available. +The **noise cutoff** metric assesses whether a unit’s spike‐amplitude distribution is truncated +at the low-end, which may be due to the high amplitude detection threshold in the deconvolution step, +i.e., if low‐amplitude spikes were missed. It does not assume a Gaussian shape; +instead, it directly compares counts in the low‐amplitude bins to counts in high‐amplitude bins. + +1. **Build a histogram** + + For each unit, divide all amplitudes into `n_bins` equally spaced bins over the range of the amplitude. + If the number of spikes is large, you may consider using a larger `n_bins`. For a small number of spikes, consider a smaller `n_bins`. + Let `n_i` denote the count in the :math:`i`-th bin. + +2. **Identify the “low” region** + - Compute the amplitude value at the specified `low_quantile` (for example, 0.10 = 10th percentile), denoted as :math:`\text{amp}_{low}`. + - Find all histogram bins whose upper edge is below that quantile value. These bins form the “low‐quantile region”. + - Compute + + .. math:: + L_{\mathrm{bins}} = \bigl\{i : \text{upper\_edge}(i) \le \text{amp}_{low}\bigr\}, \quad + \mu_{\mathrm{low}} = \frac{1}{|L_{\mathrm{bins}}|}\sum_{i\in L_{\mathrm{bins}}} n_i. + +3. **Identify the “high” region** + + - Compute the amplitude value at the specified `high_quantile` (for example, 0.25 = top 25th percentile), denoted as :math:`\text{amp}_{high}`. + - Find all histogram bins whose lower edge is greater than that quantile value. These bins form the “high‐quantile region.” + - Compute + + .. math:: + H_{\mathrm{bins}} &= \bigl\{i : \text{lower\_edge}(i)\ge \text{amp}_{high}\bigr\}, \\ + \mu_{\mathrm{high}} &= \frac{1}{|H_{\mathrm{bins}}|}\sum_{i\in H_{\mathrm{bins}}} n_i, \quad + \sigma_{\mathrm{high}} = \sqrt{\frac{1}{|H_{\mathrm{bins}}|-1}\sum_{i\in H_{\mathrm{bins}}}\bigl(n_i-\mu_{\mathrm{high}} \bigr)^2}. + +4. **Compute cutoff** + + The *cutoff* is given by how many standard deviations away the low-amplitude bins are from the high-amplitude bins, defined as + + .. math:: + \mathrm{cutoff} = \frac{\mu_{\mathrm{low}} - \mu_{\mathrm{high}}}{\sigma_{\mathrm{high}}}. + + + - If no low‐quantile bins exist, a warning is issued and `cutoff = NaN`. + - If no high‐quantile bins exist or :math:`\sigma_{\mathrm{high}} = 0`, a warning is issued and `cutoff = NaN`. + +5. **Compute the low-to-peak ratio** + + - Let :math:`M = \max_i\,n_i` be the height of the largest bin in the histogram. + - Define + + .. math:: + \mathrm{ratio} = \frac{\mu_{\mathrm{low}}}{M}. + + + - If there are no low bins, :math:`\mathrm{ratio} = NaN`. + + +Together, (cutoff, ratio) quantify how suppressed the low‐end of the amplitude distribution is relative to the top quantile and to the peak. Expectation and use ------------------- Noise cutoff attempts to describe whether an amplitude distribution is cut off. -The metric is loosely based on [Hill]_'s amplitude cutoff, but is here adapted (originally by [IBL]_) to avoid making the Gaussianity assumption on spike distributions. -Noise cutoff provides an estimate of false negative rate, so a lower value indicates fewer missed spikes (a more complete unit). +If the distribution is not truncated at the low-end, one would expect cutoff to be less than 5 and ratio less than 0.1. +It is suggested to use this metric when the amplitude histogram is **unimodal**. + +The metric is loosely based on [Hill]_'s amplitude cutoff, but is here adapted (originally by [IBL2024]_) to avoid making the Gaussianity assumption on spike distributions. + +Example code +------------ + +.. code-block:: python + + import numpy as np + import matplotlib.pyplot as plt + from spikeinterface.full as si + + # Suppose `sorting_analyzer` has been computed with spike amplitudes: + # Select units you are interested in visualizing + unit_ids = ... + + # Compute noise_cutoff: + summary_dict = compute_noise_cutoff( + sorting_analyzer=sorting_analyzer + high_quantile=0.25, + low_quantile=0.10, + n_bins=100, + unit_ids=unit_ids + ) + +Reference +--------- + +.. autofunction:: spikeinterface.qualitymetrics.misc_metrics.compute_noise_cutoff + +Examples with plots +------------------- + +Here is shown the histogram of two units, with the vertical lines separating low- and high-amplitude regions. + +- On the left, we have a unit with no truncation at the left end, and the cutoff and ratio are small. +- On the right, we have a unit with truncation at -1, and the cutoff and ratio are much larger. +.. image:: example_cutoff.png + :width: 600 Links to original implementations --------------------------------- * From `IBL implementation `_ +Note: Compared to the original implementation, we have added a comparison between the low-amplitude bins to the largest bin (`noise_ratio`). +The selection of low-amplitude bins is based on the `low_quantile` rather than the number of bins. Literature ---------- -Metric introduced by [IBL]_ (adapted from [Hill]_'s amplitude cutoff metric). +Metric introduced by [IBL2024]_. diff --git a/src/spikeinterface/qualitymetrics/misc_metrics.py b/src/spikeinterface/qualitymetrics/misc_metrics.py index a2a04a656a..95a831284e 100644 --- a/src/spikeinterface/qualitymetrics/misc_metrics.py +++ b/src/spikeinterface/qualitymetrics/misc_metrics.py @@ -37,6 +37,160 @@ _default_params = dict() +def compute_noise_cutoffs(sorting_analyzer, high_quantile=0.25, low_quantile=0.1, n_bins=100, unit_ids=None): + """ + A metric to determine if a unit's amplitude distribution is cut off as it approaches zero, without assuming a Gaussian distribution. + + Based on the histogram of the (transformed) amplitude: + + 1. This method compares counts in the lower-amplitude bins to counts in the top 'high_quantile' of the amplitude range. + It computes the mean and std of an upper quantile of the distribution, and calculates how many standard deviations away + from that mean the lower-quantile bins lie. + + 2. The method also compares the counts in the lower-amplitude bins to the count in the highest bin and return their ratio. + + Parameters + ---------- + sorting_analyzer : SortingAnalyzer + A SortingAnalyzer object. + high_quantile : float, default: 0.25 + Quantile of the amplitude range above which values are treated as "high" (e.g. 0.25 = top 25%), the reference region. + low_quantile : int, default: 0.1 + Quantile of the amplitude range below which values are treated as "low" (e.g. 0.1 = lower 10%), the test region. + n_bins: int, default: 100 + The number of bins to use to compute the amplitude histogram. + unit_ids : list or None + List of unit ids to compute the amplitude cutoffs. If None, all units are used. + + Returns + ------- + noise_cutoff_dict : dict of floats + Estimated metrics based on the amplitude distribution, for each unit ID. + + References + ---------- + Inspired by metric described in [IBL2024]_ + + """ + res = namedtuple("cutoff_metrics", ["noise_cutoff", "noise_ratio"]) + if unit_ids is None: + unit_ids = sorting_analyzer.unit_ids + + noise_cutoff_dict = {} + noise_ratio_dict = {} + if not sorting_analyzer.has_extension("spike_amplitudes"): + warnings.warn( + "`compute_noise_cutoffs` requires the 'spike_amplitudes` extension. Please run sorting_analyzer.compute('spike_amplitudes') to be able to compute `noise_cutoff`" + ) + for unit_id in unit_ids: + noise_cutoff_dict[unit_id] = np.nan + noise_ratio_dict[unit_id] = np.nan + return res(noise_cutoff_dict, noise_ratio_dict) + + amplitude_extension = sorting_analyzer.get_extension("spike_amplitudes") + peak_sign = amplitude_extension.params["peak_sign"] + if peak_sign == "both": + raise TypeError( + '`peak_sign` should either be "pos" or "neg". You can set `peak_sign` as an argument when you compute spike_amplitudes.' + ) + + amplitudes_by_units = _get_amplitudes_by_units(sorting_analyzer, unit_ids, peak_sign) + + for unit_id in unit_ids: + amplitudes = amplitudes_by_units[unit_id] + + # We assume the noise (zero values) is on the lower tail of the amplitude distribution. + # But if peak_sign == 'neg', the noise will be on the higher tail, so we flip the distribution. + if peak_sign == "neg": + amplitudes = -amplitudes + + cutoff, ratio = _noise_cutoff(amplitudes, high_quantile=high_quantile, low_quantile=low_quantile, n_bins=n_bins) + noise_cutoff_dict[unit_id] = cutoff + noise_ratio_dict[unit_id] = ratio + + return res(noise_cutoff_dict, noise_ratio_dict) + + +_default_params["noise_cutoff"] = dict(high_quantile=0.25, low_quantile=0.1, n_bins=100) + + +def _noise_cutoff(amps, high_quantile=0.25, low_quantile=0.1, n_bins=100): + """ + A metric to determine if a unit's amplitude distribution is cut off as it approaches zero, without assuming a Gaussian distribution. + + Based on the histogram of the (transformed) amplitude: + + 1. This method compares counts in the lower-amplitude bins to counts in the higher_amplitude bins. + It computes the mean and std of an upper quantile of the distribution, and calculates how many standard deviations away + from that mean the lower-quantile bins lie. + + 2. The method also compares the counts in the lower-amplitude bins to the count in the highest bin and return their ratio. + + Parameters + ---------- + amps : array-like + Spike amplitudes. + high_quantile : float, default: 0.25 + Quantile of the amplitude range above which values are treated as "high" (e.g. 0.25 = top 25%), the reference region. + low_quantile : int, default: 0.1 + Quantile of the amplitude range below which values are treated as "low" (e.g. 0.1 = lower 10%), the test region. + n_bins: int, default: 100 + The number of bins to use to compute the amplitude histogram. + + Returns + ------- + cutoff : float + (mean(lower_bins_count) - mean(high_bins_count)) / std(high_bins_count) + ratio: float + mean(lower_bins_count) / highest_bin_count + + """ + n_per_bin, bin_edges = np.histogram(amps, bins=n_bins) + + maximum_bin_height = np.max(n_per_bin) + + low_quantile_value = np.quantile(amps, q=low_quantile) + + # the indices for low-amplitude bins + low_indices = np.where(bin_edges[1:] <= low_quantile_value)[0] + + high_quantile_value = np.quantile(amps, q=1 - high_quantile) + + # the indices for high-amplitude bins + high_indices = np.where(bin_edges[:-1] >= high_quantile_value)[0] + + if len(low_indices) == 0: + warnings.warn("No bin is selected to test cutoff. Please increase low_quantile.") + return np.nan, np.nan + + # compute ratio between low-amplitude bins and the largest bin + low_counts = n_per_bin[low_indices] + mean_low_counts = np.mean(low_counts) + ratio = mean_low_counts / maximum_bin_height + + if len(high_indices) == 0: + warnings.warn("No bin is selected as the reference region. Please increase high_quantile.") + return np.nan, ratio + + if len(high_indices) == 1: + warnings.warn( + "Only one bin is selected as the reference region, and thus the standard deviation cannot be computed. " + "Please increase high_quantile." + ) + return np.nan, ratio + + # compute cutoff from low-amplitude and high-amplitude bins + high_counts = n_per_bin[high_indices] + mean_high_counts = np.mean(high_counts) + std_high_counts = np.std(high_counts) + if std_high_counts == 0: + warnings.warn("All the high-amplitude bins have the same size. Please consider changing n_bins.") + return np.nan, ratio + + cutoff = (mean_low_counts - mean_high_counts) / std_high_counts + return cutoff, ratio + + def compute_num_spikes(sorting_analyzer, unit_ids=None, **kwargs): """ Compute the number of spike across segments. @@ -805,12 +959,17 @@ def compute_amplitude_cv_metrics( def _get_amplitudes_by_units(sorting_analyzer, unit_ids, peak_sign): # used by compute_amplitude_cutoffs and compute_amplitude_medians - - if (spike_amplitudes_extension := sorting_analyzer.get_extension("spike_amplitudes")) is not None: - return spike_amplitudes_extension.get_data(outputs="by_unit", concatenated=True) + amplitudes_by_units = {} + if sorting_analyzer.has_extension("spike_amplitudes"): + spikes = sorting_analyzer.sorting.to_spike_vector() + ext = sorting_analyzer.get_extension("spike_amplitudes") + all_amplitudes = ext.get_data() + for unit_id in unit_ids: + unit_index = sorting_analyzer.sorting.id_to_index(unit_id) + spike_mask = spikes["unit_index"] == unit_index + amplitudes_by_units[unit_id] = all_amplitudes[spike_mask] elif sorting_analyzer.has_extension("waveforms"): - amplitudes_by_units = {} waveforms_ext = sorting_analyzer.get_extension("waveforms") before = waveforms_ext.nbefore extremum_channels_ids = get_template_extremum_channel(sorting_analyzer, peak_sign=peak_sign) diff --git a/src/spikeinterface/qualitymetrics/quality_metric_list.py b/src/spikeinterface/qualitymetrics/quality_metric_list.py index f7411f6376..89fea6a25e 100644 --- a/src/spikeinterface/qualitymetrics/quality_metric_list.py +++ b/src/spikeinterface/qualitymetrics/quality_metric_list.py @@ -18,6 +18,7 @@ compute_firing_ranges, compute_amplitude_cv_metrics, compute_sd_ratio, + compute_noise_cutoffs, ) from .pca_metrics import ( @@ -51,6 +52,7 @@ "firing_range": compute_firing_ranges, "drift": compute_drift_metrics, "sd_ratio": compute_sd_ratio, + "noise_cutoff": compute_noise_cutoffs, } # a dict converting the name of the metric for computation to the output of that computation @@ -81,6 +83,7 @@ "nn_noise_overlap": ["nn_noise_overlap"], "silhouette": ["silhouette"], "silhouette_full": ["silhouette_full"], + "noise_cutoff": ["noise_cutoff", "noise_ratio"], } # this dict allows us to ensure the appropriate dtype of metrics rather than allow Pandas to infer them @@ -116,4 +119,6 @@ "nn_noise_overlap": float, "silhouette": float, "silhouette_full": float, + "noise_cutoff": float, + "noise_ratio": float, } diff --git a/src/spikeinterface/qualitymetrics/tests/test_metrics_functions.py b/src/spikeinterface/qualitymetrics/tests/test_metrics_functions.py index f35ad0861c..db3e55b629 100644 --- a/src/spikeinterface/qualitymetrics/tests/test_metrics_functions.py +++ b/src/spikeinterface/qualitymetrics/tests/test_metrics_functions.py @@ -43,6 +43,7 @@ compute_quality_metrics, ) +from spikeinterface.qualitymetrics.misc_metrics import _noise_cutoff from spikeinterface.core.basesorting import minimum_spike_dtype @@ -50,6 +51,21 @@ job_kwargs = dict(n_jobs=2, progress_bar=True, chunk_duration="1s") +def test_noise_cutoff(): + """ + Generate two artifical gaussian, one truncated and one not. Check the metrics are higher for the truncated one. + """ + np.random.seed(1) + amps = np.random.normal(0, 1, 1000) + amps_trunc = amps[amps > -1] + + cutoff1, ratio1 = _noise_cutoff(amps=amps) + cutoff2, ratio2 = _noise_cutoff(amps=amps_trunc) + + assert cutoff1 <= cutoff2 + assert ratio1 <= ratio2 + + def test_compute_new_quality_metrics(small_sorting_analyzer): """ Computes quality metrics then computes a subset of quality metrics, and checks From 871c6b1191a4fe5552d8d478433241a64efdc2b0 Mon Sep 17 00:00:00 2001 From: huizizhang949 Date: Thu, 17 Jul 2025 16:21:21 +0100 Subject: [PATCH 02/18] update references --- doc/references.rst | 2 ++ 1 file changed, 2 insertions(+) diff --git a/doc/references.rst b/doc/references.rst index e4b23f1b76..16a16c1579 100644 --- a/doc/references.rst +++ b/doc/references.rst @@ -122,6 +122,8 @@ References .. [IBL] `Spike sorting pipeline for the International Brain Laboratory. 2022. `_ +.. [IBL2024] `Spike sorting pipeline for the International Brain Laboratory - Version 2. 2024. `_ + .. [Jackson] `Quantitative assessment of extracellular multichannel recording quality using measures of cluster separation. Society of Neuroscience Abstract. 2005. `_ .. [Jain] `UnitRefine: A Community Toolbox for Automated Spike Sorting Curation. 2025 `_ From 384fdc60e9443f0fc6690bb6678a774ddd1ed3e3 Mon Sep 17 00:00:00 2001 From: huizizhang949 Date: Thu, 17 Jul 2025 16:30:17 +0100 Subject: [PATCH 03/18] restore _get_amplitudes_by_units --- src/spikeinterface/qualitymetrics/misc_metrics.py | 14 ++++---------- 1 file changed, 4 insertions(+), 10 deletions(-) diff --git a/src/spikeinterface/qualitymetrics/misc_metrics.py b/src/spikeinterface/qualitymetrics/misc_metrics.py index 95a831284e..6b14975bb9 100644 --- a/src/spikeinterface/qualitymetrics/misc_metrics.py +++ b/src/spikeinterface/qualitymetrics/misc_metrics.py @@ -959,17 +959,12 @@ def compute_amplitude_cv_metrics( def _get_amplitudes_by_units(sorting_analyzer, unit_ids, peak_sign): # used by compute_amplitude_cutoffs and compute_amplitude_medians - amplitudes_by_units = {} - if sorting_analyzer.has_extension("spike_amplitudes"): - spikes = sorting_analyzer.sorting.to_spike_vector() - ext = sorting_analyzer.get_extension("spike_amplitudes") - all_amplitudes = ext.get_data() - for unit_id in unit_ids: - unit_index = sorting_analyzer.sorting.id_to_index(unit_id) - spike_mask = spikes["unit_index"] == unit_index - amplitudes_by_units[unit_id] = all_amplitudes[spike_mask] + + if (spike_amplitudes_extension := sorting_analyzer.get_extension("spike_amplitudes")) is not None: + return spike_amplitudes_extension.get_data(outputs="by_unit", concatenated=True) elif sorting_analyzer.has_extension("waveforms"): + amplitudes_by_units = {} waveforms_ext = sorting_analyzer.get_extension("waveforms") before = waveforms_ext.nbefore extremum_channels_ids = get_template_extremum_channel(sorting_analyzer, peak_sign=peak_sign) @@ -984,7 +979,6 @@ def _get_amplitudes_by_units(sorting_analyzer, unit_ids, peak_sign): return amplitudes_by_units - def compute_amplitude_cutoffs( sorting_analyzer, peak_sign="neg", From fc01b4859c0c6dde75bfdc348a1ddbd71d89c120 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Thu, 17 Jul 2025 15:30:45 +0000 Subject: [PATCH 04/18] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- src/spikeinterface/qualitymetrics/misc_metrics.py | 1 + 1 file changed, 1 insertion(+) diff --git a/src/spikeinterface/qualitymetrics/misc_metrics.py b/src/spikeinterface/qualitymetrics/misc_metrics.py index 6b14975bb9..97bf1112e7 100644 --- a/src/spikeinterface/qualitymetrics/misc_metrics.py +++ b/src/spikeinterface/qualitymetrics/misc_metrics.py @@ -979,6 +979,7 @@ def _get_amplitudes_by_units(sorting_analyzer, unit_ids, peak_sign): return amplitudes_by_units + def compute_amplitude_cutoffs( sorting_analyzer, peak_sign="neg", From 8247e413a9cb03ad898d22f0aa57c30d0d45705a Mon Sep 17 00:00:00 2001 From: huizizhang949 <79323958+huizizhang949@users.noreply.github.com> Date: Tue, 22 Jul 2025 19:18:43 +0100 Subject: [PATCH 05/18] Update doc/modules/qualitymetrics/noise_cutoff.rst Co-authored-by: Alessio Buccino --- doc/modules/qualitymetrics/noise_cutoff.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/doc/modules/qualitymetrics/noise_cutoff.rst b/doc/modules/qualitymetrics/noise_cutoff.rst index a96bd4e6c2..93d9ef29bd 100644 --- a/doc/modules/qualitymetrics/noise_cutoff.rst +++ b/doc/modules/qualitymetrics/noise_cutoff.rst @@ -5,7 +5,7 @@ Calculation ----------- -Metric describing whether an amplitude distribution is cut off as it approaches zero, similar to _amp_cutoff :ref:`amplitude cutoff ` but without a Gaussian assumption. +Metric describing whether an amplitude distribution is cut off as it approaches zero, similar to :ref:`amplitude cutoff ` but without a Gaussian assumption. The **noise cutoff** metric assesses whether a unit’s spike‐amplitude distribution is truncated at the low-end, which may be due to the high amplitude detection threshold in the deconvolution step, From d9356c89f02c6eee6dfbbc00db997336a3fc7495 Mon Sep 17 00:00:00 2001 From: huizizhang949 <79323958+huizizhang949@users.noreply.github.com> Date: Tue, 22 Jul 2025 19:18:51 +0100 Subject: [PATCH 06/18] Update doc/modules/qualitymetrics/noise_cutoff.rst Co-authored-by: Alessio Buccino --- doc/modules/qualitymetrics/noise_cutoff.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/doc/modules/qualitymetrics/noise_cutoff.rst b/doc/modules/qualitymetrics/noise_cutoff.rst index 93d9ef29bd..107f676cb3 100644 --- a/doc/modules/qualitymetrics/noise_cutoff.rst +++ b/doc/modules/qualitymetrics/noise_cutoff.rst @@ -1,4 +1,4 @@ -Noise cutoff +Noise cutoff (:code:`noise_cutoff`) ======================================== Calculation From a1890fd8fe8cc22d350dced9fa3ff528376b5a5a Mon Sep 17 00:00:00 2001 From: huizizhang949 <79323958+huizizhang949@users.noreply.github.com> Date: Tue, 22 Jul 2025 19:18:59 +0100 Subject: [PATCH 07/18] Update src/spikeinterface/qualitymetrics/misc_metrics.py Co-authored-by: Alessio Buccino --- src/spikeinterface/qualitymetrics/misc_metrics.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/spikeinterface/qualitymetrics/misc_metrics.py b/src/spikeinterface/qualitymetrics/misc_metrics.py index 97bf1112e7..4bbb202aef 100644 --- a/src/spikeinterface/qualitymetrics/misc_metrics.py +++ b/src/spikeinterface/qualitymetrics/misc_metrics.py @@ -160,7 +160,7 @@ def _noise_cutoff(amps, high_quantile=0.25, low_quantile=0.1, n_bins=100): high_indices = np.where(bin_edges[:-1] >= high_quantile_value)[0] if len(low_indices) == 0: - warnings.warn("No bin is selected to test cutoff. Please increase low_quantile.") + warnings.warn("No bin is selected to test cutoff. Please increase low_quantile. Setting noise cutoff and ratio to NaN") return np.nan, np.nan # compute ratio between low-amplitude bins and the largest bin From 62d61d8698409f162a1968e2b03b3e8d3478b322 Mon Sep 17 00:00:00 2001 From: huizizhang949 <79323958+huizizhang949@users.noreply.github.com> Date: Tue, 22 Jul 2025 19:19:05 +0100 Subject: [PATCH 08/18] Update src/spikeinterface/qualitymetrics/misc_metrics.py Co-authored-by: Alessio Buccino --- src/spikeinterface/qualitymetrics/misc_metrics.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/spikeinterface/qualitymetrics/misc_metrics.py b/src/spikeinterface/qualitymetrics/misc_metrics.py index 4bbb202aef..0b0992d7a3 100644 --- a/src/spikeinterface/qualitymetrics/misc_metrics.py +++ b/src/spikeinterface/qualitymetrics/misc_metrics.py @@ -169,7 +169,7 @@ def _noise_cutoff(amps, high_quantile=0.25, low_quantile=0.1, n_bins=100): ratio = mean_low_counts / maximum_bin_height if len(high_indices) == 0: - warnings.warn("No bin is selected as the reference region. Please increase high_quantile.") + warnings.warn("No bin is selected as the reference region. Please increase high_quantile. Setting noise cutoff to NaN") return np.nan, ratio if len(high_indices) == 1: From f78cb70013633327359754b1a1684448d8baa9a6 Mon Sep 17 00:00:00 2001 From: huizizhang949 <79323958+huizizhang949@users.noreply.github.com> Date: Tue, 22 Jul 2025 19:19:13 +0100 Subject: [PATCH 09/18] Update src/spikeinterface/qualitymetrics/misc_metrics.py Co-authored-by: Alessio Buccino --- src/spikeinterface/qualitymetrics/misc_metrics.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/spikeinterface/qualitymetrics/misc_metrics.py b/src/spikeinterface/qualitymetrics/misc_metrics.py index 0b0992d7a3..4969c205a9 100644 --- a/src/spikeinterface/qualitymetrics/misc_metrics.py +++ b/src/spikeinterface/qualitymetrics/misc_metrics.py @@ -175,7 +175,7 @@ def _noise_cutoff(amps, high_quantile=0.25, low_quantile=0.1, n_bins=100): if len(high_indices) == 1: warnings.warn( "Only one bin is selected as the reference region, and thus the standard deviation cannot be computed. " - "Please increase high_quantile." + "Please increase high_quantile. Setting noise cutoff to NaN" ) return np.nan, ratio From 6b9522c46fdd287cb78606f495e1df3d96ad9749 Mon Sep 17 00:00:00 2001 From: huizizhang949 <79323958+huizizhang949@users.noreply.github.com> Date: Tue, 22 Jul 2025 19:19:21 +0100 Subject: [PATCH 10/18] Update src/spikeinterface/qualitymetrics/misc_metrics.py Co-authored-by: Alessio Buccino --- src/spikeinterface/qualitymetrics/misc_metrics.py | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/src/spikeinterface/qualitymetrics/misc_metrics.py b/src/spikeinterface/qualitymetrics/misc_metrics.py index 4969c205a9..b587a5f7d9 100644 --- a/src/spikeinterface/qualitymetrics/misc_metrics.py +++ b/src/spikeinterface/qualitymetrics/misc_metrics.py @@ -184,7 +184,10 @@ def _noise_cutoff(amps, high_quantile=0.25, low_quantile=0.1, n_bins=100): mean_high_counts = np.mean(high_counts) std_high_counts = np.std(high_counts) if std_high_counts == 0: - warnings.warn("All the high-amplitude bins have the same size. Please consider changing n_bins.") + warnings.warn( + "All the high-amplitude bins have the same size. Please consider changing n_bins. " + "Setting noise cutoff to NaN" + ) return np.nan, ratio cutoff = (mean_low_counts - mean_high_counts) / std_high_counts From be50e145b48ddc4b8786f8ffe0386491d8567737 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Tue, 22 Jul 2025 18:20:27 +0000 Subject: [PATCH 11/18] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- src/spikeinterface/qualitymetrics/misc_metrics.py | 8 ++++++-- 1 file changed, 6 insertions(+), 2 deletions(-) diff --git a/src/spikeinterface/qualitymetrics/misc_metrics.py b/src/spikeinterface/qualitymetrics/misc_metrics.py index b587a5f7d9..e2102d3ff2 100644 --- a/src/spikeinterface/qualitymetrics/misc_metrics.py +++ b/src/spikeinterface/qualitymetrics/misc_metrics.py @@ -160,7 +160,9 @@ def _noise_cutoff(amps, high_quantile=0.25, low_quantile=0.1, n_bins=100): high_indices = np.where(bin_edges[:-1] >= high_quantile_value)[0] if len(low_indices) == 0: - warnings.warn("No bin is selected to test cutoff. Please increase low_quantile. Setting noise cutoff and ratio to NaN") + warnings.warn( + "No bin is selected to test cutoff. Please increase low_quantile. Setting noise cutoff and ratio to NaN" + ) return np.nan, np.nan # compute ratio between low-amplitude bins and the largest bin @@ -169,7 +171,9 @@ def _noise_cutoff(amps, high_quantile=0.25, low_quantile=0.1, n_bins=100): ratio = mean_low_counts / maximum_bin_height if len(high_indices) == 0: - warnings.warn("No bin is selected as the reference region. Please increase high_quantile. Setting noise cutoff to NaN") + warnings.warn( + "No bin is selected as the reference region. Please increase high_quantile. Setting noise cutoff to NaN" + ) return np.nan, ratio if len(high_indices) == 1: From bfb5fe70d373ec04a6d65961651d063414db3ebc Mon Sep 17 00:00:00 2001 From: huizizhang949 Date: Tue, 22 Jul 2025 19:49:40 +0100 Subject: [PATCH 12/18] update doc --- doc/modules/qualitymetrics/noise_cutoff.rst | 13 ++++++++----- 1 file changed, 8 insertions(+), 5 deletions(-) diff --git a/doc/modules/qualitymetrics/noise_cutoff.rst b/doc/modules/qualitymetrics/noise_cutoff.rst index 107f676cb3..1d2a1410ed 100644 --- a/doc/modules/qualitymetrics/noise_cutoff.rst +++ b/doc/modules/qualitymetrics/noise_cutoff.rst @@ -1,11 +1,11 @@ -Noise cutoff (:code:`noise_cutoff`) +Noise cutoff ======================================== Calculation ----------- -Metric describing whether an amplitude distribution is cut off as it approaches zero, similar to :ref:`amplitude cutoff ` but without a Gaussian assumption. +Metric describing whether an amplitude distribution is cut off as it approaches zero, similar to _amp_cutoff :ref:`amplitude cutoff ` but without a Gaussian assumption. The **noise cutoff** metric assesses whether a unit’s spike‐amplitude distribution is truncated at the low-end, which may be due to the high amplitude detection threshold in the deconvolution step, @@ -24,7 +24,7 @@ instead, it directly compares counts in the low‐amplitude bins to counts in hi - Compute .. math:: - L_{\mathrm{bins}} = \bigl\{i : \text{upper\_edge}(i) \le \text{amp}_{low}\bigr\}, \quad + L_{\mathrm{bins}} = \bigl\{i : \text{upper_edge}_i \le \text{amp}_{low}\bigr\}, \quad \mu_{\mathrm{low}} = \frac{1}{|L_{\mathrm{bins}}|}\sum_{i\in L_{\mathrm{bins}}} n_i. 3. **Identify the “high” region** @@ -34,7 +34,7 @@ instead, it directly compares counts in the low‐amplitude bins to counts in hi - Compute .. math:: - H_{\mathrm{bins}} &= \bigl\{i : \text{lower\_edge}(i)\ge \text{amp}_{high}\bigr\}, \\ + H_{\mathrm{bins}} &= \bigl\{i : \text{lower_edge}_i \ge \text{amp}_{high}\bigr\}, \\ \mu_{\mathrm{high}} &= \frac{1}{|H_{\mathrm{bins}}|}\sum_{i\in H_{\mathrm{bins}}} n_i, \quad \sigma_{\mathrm{high}} = \sqrt{\frac{1}{|H_{\mathrm{bins}}|-1}\sum_{i\in H_{\mathrm{bins}}}\bigl(n_i-\mu_{\mathrm{high}} \bigr)^2}. @@ -67,7 +67,10 @@ Expectation and use ------------------- Noise cutoff attempts to describe whether an amplitude distribution is cut off. -If the distribution is not truncated at the low-end, one would expect cutoff to be less than 5 and ratio less than 0.1. +Larger values of `cutoff` and `ratio` suggest that the distribution is cut-off. +IBL uses the default value of 1 to choose the number of lower bins, with a suggested threshold of 5 for `cutoff` and 0.1 for `ratio` to determine whether a unit is cut off or not. +In practice, the IBL threshold is quite conservative, and a lower threshold might be better for your data. +We suggest plotting the data using the `plot_amplitudes` widget to view your data when choosing your threshold. It is suggested to use this metric when the amplitude histogram is **unimodal**. The metric is loosely based on [Hill]_'s amplitude cutoff, but is here adapted (originally by [IBL2024]_) to avoid making the Gaussianity assumption on spike distributions. From a44c0fafe7548171549baf85551c661bbdb23433 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Tue, 22 Jul 2025 18:50:17 +0000 Subject: [PATCH 13/18] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- doc/modules/qualitymetrics/noise_cutoff.rst | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/doc/modules/qualitymetrics/noise_cutoff.rst b/doc/modules/qualitymetrics/noise_cutoff.rst index 1d2a1410ed..e65e80dae6 100644 --- a/doc/modules/qualitymetrics/noise_cutoff.rst +++ b/doc/modules/qualitymetrics/noise_cutoff.rst @@ -67,9 +67,9 @@ Expectation and use ------------------- Noise cutoff attempts to describe whether an amplitude distribution is cut off. -Larger values of `cutoff` and `ratio` suggest that the distribution is cut-off. -IBL uses the default value of 1 to choose the number of lower bins, with a suggested threshold of 5 for `cutoff` and 0.1 for `ratio` to determine whether a unit is cut off or not. -In practice, the IBL threshold is quite conservative, and a lower threshold might be better for your data. +Larger values of `cutoff` and `ratio` suggest that the distribution is cut-off. +IBL uses the default value of 1 to choose the number of lower bins, with a suggested threshold of 5 for `cutoff` and 0.1 for `ratio` to determine whether a unit is cut off or not. +In practice, the IBL threshold is quite conservative, and a lower threshold might be better for your data. We suggest plotting the data using the `plot_amplitudes` widget to view your data when choosing your threshold. It is suggested to use this metric when the amplitude histogram is **unimodal**. From 36b19012fd841ef141bf1eaa7d2e8525030ddfca Mon Sep 17 00:00:00 2001 From: huizizhang949 Date: Tue, 22 Jul 2025 20:49:11 +0100 Subject: [PATCH 14/18] update latex --- doc/modules/qualitymetrics/noise_cutoff.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/doc/modules/qualitymetrics/noise_cutoff.rst b/doc/modules/qualitymetrics/noise_cutoff.rst index e65e80dae6..fccc43f5b6 100644 --- a/doc/modules/qualitymetrics/noise_cutoff.rst +++ b/doc/modules/qualitymetrics/noise_cutoff.rst @@ -16,7 +16,7 @@ instead, it directly compares counts in the low‐amplitude bins to counts in hi For each unit, divide all amplitudes into `n_bins` equally spaced bins over the range of the amplitude. If the number of spikes is large, you may consider using a larger `n_bins`. For a small number of spikes, consider a smaller `n_bins`. - Let `n_i` denote the count in the :math:`i`-th bin. + Let :math:`n_i` denote the count in the :math:`i`-th bin. 2. **Identify the “low” region** - Compute the amplitude value at the specified `low_quantile` (for example, 0.10 = 10th percentile), denoted as :math:`\text{amp}_{low}`. From 1196efda41508a5f207bd829f3168f2f956d6165 Mon Sep 17 00:00:00 2001 From: huizizhang949 Date: Tue, 22 Jul 2025 20:53:29 +0100 Subject: [PATCH 15/18] update doc --- doc/modules/qualitymetrics/noise_cutoff.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/doc/modules/qualitymetrics/noise_cutoff.rst b/doc/modules/qualitymetrics/noise_cutoff.rst index fccc43f5b6..d7713183fb 100644 --- a/doc/modules/qualitymetrics/noise_cutoff.rst +++ b/doc/modules/qualitymetrics/noise_cutoff.rst @@ -5,7 +5,7 @@ Calculation ----------- -Metric describing whether an amplitude distribution is cut off as it approaches zero, similar to _amp_cutoff :ref:`amplitude cutoff ` but without a Gaussian assumption. +Metric describing whether an amplitude distribution is cut off as it approaches zero, similar to :ref:`amplitude cutoff ` but without a Gaussian assumption. The **noise cutoff** metric assesses whether a unit’s spike‐amplitude distribution is truncated at the low-end, which may be due to the high amplitude detection threshold in the deconvolution step, From 5803fdb26b693e30e0b39c9f112450af98b512f3 Mon Sep 17 00:00:00 2001 From: huizizhang949 Date: Tue, 22 Jul 2025 20:55:15 +0100 Subject: [PATCH 16/18] doc --- doc/modules/qualitymetrics/noise_cutoff.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/doc/modules/qualitymetrics/noise_cutoff.rst b/doc/modules/qualitymetrics/noise_cutoff.rst index d7713183fb..673f03ebb1 100644 --- a/doc/modules/qualitymetrics/noise_cutoff.rst +++ b/doc/modules/qualitymetrics/noise_cutoff.rst @@ -1,4 +1,4 @@ -Noise cutoff +Noise cutoff (:code:`noise_cutoff`) ======================================== Calculation From d6b971f8636777f5445fd5b37cf5cfa974f917d2 Mon Sep 17 00:00:00 2001 From: Alessio Buccino Date: Wed, 23 Jul 2025 17:41:18 +0200 Subject: [PATCH 17/18] Update doc/modules/qualitymetrics/noise_cutoff.rst --- doc/modules/qualitymetrics/noise_cutoff.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/doc/modules/qualitymetrics/noise_cutoff.rst b/doc/modules/qualitymetrics/noise_cutoff.rst index 673f03ebb1..44a5ba1741 100644 --- a/doc/modules/qualitymetrics/noise_cutoff.rst +++ b/doc/modules/qualitymetrics/noise_cutoff.rst @@ -1,5 +1,5 @@ Noise cutoff (:code:`noise_cutoff`) -======================================== +=================================== Calculation ----------- From 1a4a1dd06cd8f0d5abfb0e08c4816674a0673070 Mon Sep 17 00:00:00 2001 From: chrishalcrow Date: Thu, 24 Jul 2025 10:50:22 +0100 Subject: [PATCH 18/18] docs formatting --- doc/modules/qualitymetrics/noise_cutoff.rst | 33 +++++++++++---------- 1 file changed, 17 insertions(+), 16 deletions(-) diff --git a/doc/modules/qualitymetrics/noise_cutoff.rst b/doc/modules/qualitymetrics/noise_cutoff.rst index 44a5ba1741..af15571cf2 100644 --- a/doc/modules/qualitymetrics/noise_cutoff.rst +++ b/doc/modules/qualitymetrics/noise_cutoff.rst @@ -14,13 +14,13 @@ instead, it directly compares counts in the low‐amplitude bins to counts in hi 1. **Build a histogram** - For each unit, divide all amplitudes into `n_bins` equally spaced bins over the range of the amplitude. - If the number of spikes is large, you may consider using a larger `n_bins`. For a small number of spikes, consider a smaller `n_bins`. + For each unit, divide all amplitudes into ``n_bins`` equally spaced bins over the range of the amplitude. + If the number of spikes is large, you may consider using a larger ``n_bins``. For a small number of spikes, consider a smaller ``n_bins``. Let :math:`n_i` denote the count in the :math:`i`-th bin. 2. **Identify the “low” region** - - Compute the amplitude value at the specified `low_quantile` (for example, 0.10 = 10th percentile), denoted as :math:`\text{amp}_{low}`. - - Find all histogram bins whose upper edge is below that quantile value. These bins form the “low‐quantile region”. + - Compute the amplitude value at the specified ``low_quantile`` (for example, 0.10 = 10th percentile), denoted as :math:`\text{amp}_{low}`. + - Find all histogram bins whose upper edge is below that quantile value. These bins form the "low‐quantile region". - Compute .. math:: @@ -29,8 +29,8 @@ instead, it directly compares counts in the low‐amplitude bins to counts in hi 3. **Identify the “high” region** - - Compute the amplitude value at the specified `high_quantile` (for example, 0.25 = top 25th percentile), denoted as :math:`\text{amp}_{high}`. - - Find all histogram bins whose lower edge is greater than that quantile value. These bins form the “high‐quantile region.” + - Compute the amplitude value at the specified ``high_quantile`` (for example, 0.25 = top 25th percentile), denoted as :math:`\text{amp}_{high}`. + - Find all histogram bins whose lower edge is greater than that quantile value. These bins form the "high‐quantile region". - Compute .. math:: @@ -46,8 +46,8 @@ instead, it directly compares counts in the low‐amplitude bins to counts in hi \mathrm{cutoff} = \frac{\mu_{\mathrm{low}} - \mu_{\mathrm{high}}}{\sigma_{\mathrm{high}}}. - - If no low‐quantile bins exist, a warning is issued and `cutoff = NaN`. - - If no high‐quantile bins exist or :math:`\sigma_{\mathrm{high}} = 0`, a warning is issued and `cutoff = NaN`. + - If no low‐quantile bins exist, a warning is issued and ``cutoff = NaN``. + - If no high‐quantile bins exist or :math:`\sigma_{\mathrm{high}} = 0`, a warning is issued and ``cutoff = NaN``. 5. **Compute the low-to-peak ratio** @@ -61,19 +61,20 @@ instead, it directly compares counts in the low‐amplitude bins to counts in hi - If there are no low bins, :math:`\mathrm{ratio} = NaN`. -Together, (cutoff, ratio) quantify how suppressed the low‐end of the amplitude distribution is relative to the top quantile and to the peak. +Together, ``(cutoff, ratio)`` quantify how suppressed the low‐end of the amplitude distribution is relative to the top quantile and to the peak. Expectation and use ------------------- Noise cutoff attempts to describe whether an amplitude distribution is cut off. -Larger values of `cutoff` and `ratio` suggest that the distribution is cut-off. -IBL uses the default value of 1 to choose the number of lower bins, with a suggested threshold of 5 for `cutoff` and 0.1 for `ratio` to determine whether a unit is cut off or not. -In practice, the IBL threshold is quite conservative, and a lower threshold might be better for your data. -We suggest plotting the data using the `plot_amplitudes` widget to view your data when choosing your threshold. +Larger values of ``cutoff`` and ``ratio`` suggest that the distribution is cut-off. +IBL uses the default value of 1 (equivalent to e.g. ``low_quantile=0.01, n_bins=100``) to choose the number of +lower bins, with a suggested threshold of 5 for ``cutoff`` to determine whether a unit is cut off or not. +In practice, the IBL threshold is quite conservative, and a lower threshold might work better for your data. +We suggest plotting the data using the :py:func:`~spikeinterface.widgets.plot_amplitudes` widget to view your data when choosing your threshold. It is suggested to use this metric when the amplitude histogram is **unimodal**. -The metric is loosely based on [Hill]_'s amplitude cutoff, but is here adapted (originally by [IBL2024]_) to avoid making the Gaussianity assumption on spike distributions. +The metric is loosely based on [Hill]_'s amplitude cutoff, but is here adapted (originally by [IBL2024]_) to avoid making the Gaussian assumption on spike distributions. Example code ------------ @@ -118,8 +119,8 @@ Links to original implementations * From `IBL implementation `_ -Note: Compared to the original implementation, we have added a comparison between the low-amplitude bins to the largest bin (`noise_ratio`). -The selection of low-amplitude bins is based on the `low_quantile` rather than the number of bins. +Note: Compared to the original implementation, we have added a comparison between the low-amplitude bins to the largest bin (``noise_ratio``). +The selection of low-amplitude bins is based on the ``low_quantile`` rather than the number of bins. Literature ----------