PRX Life 進(jìn)展：細(xì)胞狀態(tài)轉(zhuǎn)化中基因調(diào)控網(wǎng)絡(luò)的動(dòng)力學(xué)研究 Transiently Increased Coordination in Gene Regulation During Cell Phenotypic Transitions

PRXLIFE2,043009(2024)

TransientlyIncreasedCoordinationinGeneRegulationDuringCellPhenotypicTransitions

WeikangWang

,

1,2,*

KeNi

,

3

,4

DantePoe

,

3

,4

andJianhuaXing4,

5,

6

,?

1CASKeyLaboratoryofTheoreticalPhysics,

InstituteofTheoreticalPhysics,

ChineseAcademyofSciences,

Beijing100190,China

2SchoolofPhysicalSciences,

UniversityofChineseAcademyofSciences,

Beijing100049,China

3JointCarnegieMellonUniversity-

UniversityofPittsburgh

Ph.D.PrograminComputationalBiology,Pittsburgh,Pennsylvania15213,USA

4DepartmentofComputationalandSystemsBiology,SchoolofMedicine,

UniversityofPittsburgh,

Pittsburgh,Pennsylvania15213,USA

5DepartmentofPhysicsandAstronomy,

UniversityofPittsburgh,

Pittsburgh,Pennsylvania15232,USA

6UPMC–HillmanCancerCenter,

UniversityofPittsburgh,

Pittsburgh,Pennsylvania15232,USA

畫

(Received11January2024;accepted30September2024;published5November2024)

Phenotypetransitionsoccurinmanybiologicalprocessessuchasdifferentiationandreprogramming.Afun-damentalquestionishowcellscoordinateswitchingofgeneexpressionclusters.Byanalyzingsingle-cellRNAsequencingdatawithintheframeworkoftransitionpaththeory,westudiedthegenome-wideexpressionprogramswitchingin?vedifferentcelltransitionprocesses.Foreachprocesswereconstructedareactioncoordinatedescribingthetransitionprogression,andweinferredthegeneregulatorynetworkalongthisreactioncoordinate.Inallprocessesweobservedacommonpattern:theoveralleffectivenumberandstrengthofregulationbetweendifferentcommunitiesincrease?rstandthendecrease.Thischangeisaccompaniedbysimilarchangesingeneregulatorynetworkfrustration—de?nedastheoverallcon?ictbetweentheregulationreceivedbygenesandtheirexpressionstates.Complementingpreviousstudiessuggestingthatbiologicalnetworksaremodularizedtocontainperturbationeffectslocally,ouranalysesonthe?vecelltransitionprocesseslikelyrevealageneralprinciple:duringacellphenotypictransition,intercommunityinteractionsincreasetoconcertedlycoordinateglobalgeneexpressionreprogrammingandcanalizetospeci?ccellphenotype,asWaddingtonvisioned.

DOI:

10.1103/PRXLife.2.043009

I.INTRODUCTION

Alastingtopicinscienceandengineeringisunderstandinghowadynamicalsystemtransitsfromonestableattractortoanewoneinthecorrespondingstatespace[

1]

.Oneimportanttypeoftransitionthathasgarneredincreasedinterestrecentlyisthetransitionbetweendifferentcellphenotypes.Thisin-terestispartlydrivenbytheavailabilityofgenome-widecharacterizationsofcellgeneexpressionstatesthroughoutthetransitionprocess,facilitatedbyadvancesinsingle-cellgenomicstechniques[

2

–

4]

.Acellisanonlineardynamicalsystemgovernedbyacomplexregulatorynetworkformedbymanyinteractinggenes,whichcanhavemultiplestableattractorscorrespondingtodifferentcellphenotypes

[5,6]

.Typically,alargenumberofphenotype-speci?cgenesmain-tainaspeci?cphenotypethroughmutualactivationwhilesuppressingexpressionofgenescorrespondingtootherexclu-sivephenotypes[

7]

.Insomesenseitresemblesaspinsystemsegregatingintoupwardanddownwarddomains.Whenacellphenotypictransitionoccurs,thegenesneedtoswitchtheir

*Contactauthor:wangwk@?Contactauthor:xing1@

PublishedbytheAmericanPhysicalSocietyunderthetermsofthe

CreativeCommonsAttribution4.0International

license.Furtherdistributionofthisworkmustmaintainattributiontotheauthor(s)andthepublishedarticle’stitle,journalcitation,andDOI.

expressionstatus,analogousto?ippingsomeupwardanddownwardspindomains[

8

,9]

.

Akeyquestionishowacellphenotypictransition,oracellstatetransitionprocessingeneral,proceeds.Answeringthisquestionrequiresexamininghowagenome-widegeneregulatorynetworkchangesduringacellphenotypictransi-tion[Fig.

1(a)]

.Thetransitionmaybesequentialinvolvingthedeactivationoflinks(e.g.,betweenpairsofgenes)intheregulatorynetworkthatde?nestheinitialphenotype,followedbytheactivationoflinksintheregulatorynetworkthatde?nesthe?nalphenotype.Forexample,ithasbeenobservedthatcellfatedecision-makingpreferentiallyoccursaftercellspassthroughtheMphase,wheretranscriptionfactorsdisassociatefromthecondensedchromatin,allowingtheresetofthecell-type-speci?cexpressionprogramatthetranscriptionlevel;then,intheG1phase,cellsrecruittranscriptionfactorstoactivatetheexpressionofgenesofthenewcelltype[

10]

.Alternatively,inaconcertedtransition,thedeactivationandactivationofregulatorynetworklinksmayhappenconcur-

rently[4,7,11]

.Figure

1(a)

summarizesqualitativelydifferentcharacteristicsofthesetwomechanisms.

Developmentsofsingle-cellRNAsequencing(scRNA-seq),alongwithadvancedanalysistools,havesigni?cantlyexpandedourunderstandingofcellphenotypictransitions

[12–

15]

.Duetoitsdestructivenature,scRNA-seqdataonlyprovidesnapshotinformationofcellstates.Consequently,pseudotimeanalysishasbeenwidelyusedtoorderthesam-pledcellstatesbasedontheirexpressionpro?les,utilizingtoolssuchasMonocle,Scanpy,andSeurat[

2

,12,14,16–

21]

.ArecentlydevelopedRNAvelocityformalism[

22

]leverages

2835-8279/2024/2(4)/043009(15)043009-1PublishedbytheAmericanPhysicalSociety

WANG,NI,POE,ANDXINGPRXLIFE2,043009(2024)

043009-2

FIG.1.Overviewoftheanalysispipeline.(a)Schematicplotofdynamicsofcellstateandunderlyinggeneregulatorynetworksduringcellphenotypictransition.Filledcirclesandemptycirclerepresentactivegenesandsilentgenesseparately.Colorsindicatemarkergenesofdifferentcellstates.Arrowsrepresentactivation,whilebluntarrowsrepresentinhibition.(b)Flowchartofsingle-celltrajectorysimulationandreactioncoordinatewithRNAvelocityanalysis.Left:single-cellRNAsequencingdata;middle:thecoloreddotswithanarrowarecellswiththeirRNAvelocities(differentcolorsindicatedifferentcelltypes);right:single-celltrajectorysimulation(cyanline)andreactioncoordinate(largecoloreddots)ofthecellphenotypictransitionprocess.GraydotsrepresentsinglecellandwhitelinesareboundariesofVoronoigridsgeneratedbyreactioncoordinatepoints.(c)Schematicplotofgeneregulatorynetworkanalyses.ThegeneregulatorynetworkisinferredfromRNAvelocityofscRNA-seqdata.Binarizationofgenesallowsanalysesonthevariationofgeneregulatorynetworkcommunity(left)andfrustration(right)duringcellphenotypictransition.Left:communityanalysisofgeneregulatorynetwork.Thedashedellipsesrepresentcommunities.Colorsindicatemarkergenesofdifferentcellstates(light-bluecirclesindicategenesareactiveininitialstates,andpurplecirclesindicategenesareactivein?nalstate).Right:frustrationingeneregulatorynetwork.Filledcircles(upwardsarrow)representactivegenes.Emptycircle(downwardsarrow)representsasilentgene.

thecopynumbersofunsplicedandsplicedmessengerRNAwithinasinglecelltoinferthetimederivativeofacellstateinthegeneexpressionspace,predictingthesingle-cellstateoverthenextfewhours.RNAvelocityhasbeenemployedforcellstatepredictionandhasrevealeddynamicinformationaboutcellphenotypictransitions,suchasdifferentiation.

TheseapplicationsofRNAvelocitieshaveinspirednumer-ousdevelopmentsaimedatfurtherimprovingtheaccuracyofRNAvelocityestimation[

23

–

26]

.

UsingvariousformsofRNAvelocitydata,Qiuetal.developedageneralDynamoframeworkthatintegratessingle-celldataanalyseswithcontinuousdynamicalsystems

TRANSIENTLYINCREASEDCOORDINATIONINGENE…PRXLIFE2,043009(2024)

043009-3

functionF(x)fromsnapshotsingle-celldataandperform

theory-basedmechanisticmodelingofcellularsystems[

27]

.Zhangetal.furtherdevelopedacomplementarydiscretemodelingframework,Graph-Dynamo,whichisconvenientfornumericalanalyses[

28]

.Inthiswork,byapplyingthesemodelingframeworksalongwithreactionratetheoriesandnetworksciencetheories,weanalyzedscRNA-seqdatasetsofseveralcellphenotypictransitionprocesses,andweuncov-eredsharedprinciplesonhowtheglobalreprogrammingofgeneexpressionproceeds.

II.RESULTS

A.AnalyzescRNA-seqdatawithintheframework

ofdynamicalsystemstheory

Dynamoisageneralframeworkforreconstructingaset

ofdynamicalequationsthatdescribehowacellstateevolves

overtimeusingsingle-celldataandtheinstantaneousve-

locitiesofcellstatechanges[

27]

.Inthisstudy,weused

splicing-basedRNAvelocities[

22

][seeMaterialsandMeth-

odsintheSupplementalMaterial(SM)[

29

]].Denotexas

acellstatevector,andassumethatitstemporalevolution

dt=F(x).WithDynamo,onecanconstructthevector?eld

xdescribedbyasetofdynamicalequationsgenericallyas

trajectorysimulationsusingthisfunction.

Topropagatesingle-celltrajectoriesinthehigh-dimensionalstatespace,weadoptedGraph-Dynamo[

28]

inthisstudytodiscretizethecontinuousdynamicalequationsasMarkoviantransitionsonadiscretegraphformedbysingle-celldatapoints(MaterialsandMethodsintheSM

[29])

.ThismethodutilizestheFokker-Planckequationformalismtomodelcellstatetransitionsonthedata-formedgraphnetwork,whichisinvariantunderrepresentationtransformationandpreservesthetopologicalanddynamicalpropertiesofthesystemdynamics.Thisdiscretizationcanberegardedasanextensionofdiscretizingdynamicalequationsonaregularlattice[

30

,31],typicallyusedinnumerical

simulations,toirregularlattices.Thelatterapproachleveragesthefactthatsingle-celldataresideinalow-dimensionalmanifoldthatcanbediscretizedbythedatapoints,bypassingthecurse-of-dimensionalityifaregularlatticeisused.Wethenappliedageneralized?nite-temperaturestringmethod(MaterialsandMethodsintheSM[

29

])[

32

]toanalyzethesimulatedsingle-celltrajectories.Thestringmethod,widelyusedforanalyzingtransitionpathsinchemicalreactions,essentiallyprovidesareactioncoordinate,acentralconceptinratetheories

[1]thatusesaone-dimensionalmanifold

tore?ectprogressionofthereactionprocess.Noticingthesimilaritybetweencellphenotypetransitionsandchemicalreactions[

7

],weutilizedareactioncoordinatetostudycellphenotypictransitions(SM2[

29

]).Weidenti?edaone-dimensionalreactioncoordinate[denotedby[r],whereristheindexofthediscretizedreactioncoordinate(usingaVoronoigrid)]tocharacterizetheprogressionofacellphenotypictransitionprocess[Fig.

1(b)]

.Inthecontextofcellphenotypictransitions,areactioncoordinatecanbeconceptuallyrelatedtoapseudotimetrajectoryusedinthesingle-cellgenomics?eld[

12

,14]

.

Withthereactioncoordinateidenti?ed,weaimedtoun-ravelhowthegeneexpressionprogramchangesalongthereactioncoordinateduringacellphenotypictransitionpro-cess.Toachievethisaim,wedevelopedaproceduretocoarse-grainacontinuousvector?eldintoadirectedgenereg-ulatorynetwork.First,tofacilitateanalysisoftheregulationprogramalongthereactioncoordinate,weemployedapartialleast-squaresregression(PLSR)combinedwithalocalfalsediscoverrate(LFDR)methodtocoarse-grainagenerallynon-linearvectorfunctionintoalinearmodel,Fi(x)→ΣjFijxj(MaterialsandMethodsintheSM[

29

]).HereFisamatrixwithFij>0ifgenejactivatesgenei,0fornoregulation,and<0forinhibition.Next,weidenti?edgenesshowingaswitchlikeexpressionchangeduringthephenotypetransition.Wefurthercoarse-grainedthecontinuousgeneexpressionstateofsuchageneiintodiscretestates,si=0forsilenceand1foractiveexpression[Fig.

1(c)]

.

Followingtheaboveprocedure,weobtainedadirectedgeneregulatorynetworkmodelforagivensystem,whichallowedapplicationofvariousmethodsdevelopedinnetworkscienceanalyses,suchascommunitydetection[Fig.

1(c)]

(MaterialsandMethodsintheSM[

29

]).Speci?cally,wede?nedthataneffectiveedgeexistsbetweenaregulatorygenejandatargetgeneiatacellstateifandonlyifFij0andsj=1.Thatis,theeffectivegeneregulatorynetworkiscell-state-speci?c.Wealsointroducedaconceptthataregulationisfrustratedatacellstateiftheregulationofgenej(withitsownexpressionstatesj=1)ongeneicontradictstheexpressionstateofgenei,i.e.,Fij<0andsi=1,orFij>0andsi=0.

B.DynamicalmodelreconstructedfromscRNA-seqdatadescribesdevelopmentofgranulelineageindentategyrus

Thedentategyrusislocatedinthehippocampusofthebrainandiscrucialformemoryformation[

33

,34]

.Duringneurogenesis,radialglia-likecellsdifferentiatethroughneuralintermediateprogenitorcells(nIPCs),neuroblast1and2,immaturegranulecells,andeventuallyintomaturegranulecells[Fig.

2(a)

][35]

.Severaltranscriptionfactors,suchasHes5,Sox2,andCREB,playimportantrolesatvariousstepsofthisprocess[

36]

.

UsingthedatasetofHochgerneretal.

[35]asinput,

weobtainedanarrayofreactioncoordinatepointsfromthesimulatedtrajectories[Fig.

2(b)

],whichre?ectthecontin-uousprogressionfromradialgliatogranulecellsthroughasetofintermediatecelltypes[Fig.

2(c)

][32,37,38]

.Fol-lowingthecoarse-grainingprocedure,weobtainedregulatoryrelationshipsamonggenes.TableS1listsexistinglitera-turesupportforsomeoftheinferredregulations,validatingourprocedureofinferringregulationrelationshipsfromthescRNA-seqdata.Forinstance,weidenti?edtheregulationofSMAD1onCD9,whichisakeypathwayinthedifferentiationprocess[

39]

.

C.Networkanalysesrevealincreaseoffrustrationandnetwork

heterogeneityduringtransition

Aconcerted,ratherthansequential,mechanismpredictsatransientincreaseingene-geneinteractions[Fig.

1(a)]

.Fordentategyrusdevelopment,thenumberofeffectiveedges

WANG,NI,POE,ANDXINGPRXLIFE2,043009(2024)

043009-4

FIG.2.AnalysesofscRNA-seqdataofdentategyrusneurogenesisrevealaconcertedtransitionmechanism.(a)scRNA-seqdataandRNAvelocity-basedtransitiongraphshowninthecellexpressionstatespace[showninthe2Dleadingprincipalcomponents(PCs)space].Eachdotrepresentsacell,andeachedgebetweentwodotsindicatesatransitionbetweencellstatescorrespondingtothetwocells.Colorrepresentscelltype.Arrowindicatesdirectionoftransition.(b)Atypicalsingle-celltrajectorysimulatedonthetransitiongraph,illustratedinthe2DleadingPCsspace.ThetrajectorystartsfromtheRadial-Gliacelltypeandtransitsintothegranulecelltype.Dotcolorrepresentstheprogressionofthetrajectory(frombluetored).(c)1Dreactioncoordinatereconstructedfromthesimulatedsingle-celltrajectories.Largecoloreddotsrepresentthereactioncoordinatepoints(startsfromblueandendsinred).Voronoigridsaregeneratedwiththereactioncoordinatepoints.Thesmalldotsarecellswithcolorindicatingcelltype.(d)Frustrationscorealongthereactioncoordinateofdentategyrusneurogenesis.Themeanandvarianceateachreactioncoordinatepointwerecalculatedusingk-nearest-neighboringcellsofthecorrespondingreactioncoordinatepointwithintheVoronoigrid.Thegreenviolinrepresentsdistributionandthedashedlineindicatestheaveragefrustrationscoreofrandomsamples.Notethatsomerandomstatesmaynotbebiologicallyaccessible.(e)Variationofdistributionofstateoverlapalongthereactioncoordinateofdentategyrusneurogenesis.Colorsofthedistributionscorrespondtothatofreactioncoordinatepoints.Thedistributionsofinitialand?nalreactioncoordinateareplottedwithadottedline.InsertedgraphistheFWHMofthedistributionsalongthereactioncoordinate.(f)Cell-speci?cvariationofthenumberofeffectiveintercommunityedges(representedwithcolor)inthegeneregulatorynetworkofeachsinglecell.(g)Evolutionofthenumberofeffectiveintercommunityedgesalongthereactioncoordinateduringdentategyrusneurogenesis.Eachnoderepresentsagene.Thecolorofthenodeindicatestheindexofthecommunity.Arrowsrepresentdirectionofregulation.ristheindexinthereactioncoordinate.(h)Schematicoftheconcertedmechanismforacellphenotypictransition.Thedashed-lineboxesrepresentcommunities.Filledcirclesrepresentactivegenesandemptycirclesrepresentsilentgenes.

TRANSIENTLYINCREASEDCOORDINATIONINGENE…PRXLIFE2,043009(2024)

043009-5

indeedincreasedinitiallyandthendecreasedalongthereactioncoordinate(Fig.S1a),supportingtheconcertedmechanismforthisprocess.

Inaconcertedmechanism,oneexpectsatransientcoex-pressionofgenesthatdonotcoexpressinastablephenotype,leadingtocon?ictbetweentheexpressionstateofageneandtheregulationactingonit.Weusedfrustrationtoquantifysuchcon?ictingregulations.Wede?nedtheoverallfrustrationscoreofacell-speci?cgeneregulatorynetworkasthefractionoffrustratededgesoutofalledgesintheentirenetworkofthecell[

40]

.Fordevelopmentofthegranulelineage,theaveragefrustrationscorealongthereactioncoordinateincreases?rst,reachingapeakcorrespondingtoneuroblastcells,andthendecreases[Fig.

2(d)

],consistentwiththeconcertedmecha-nismratherthanthesequentialmechanism.Forcomparison,wealsocalculatedthefrustrationscoreofrandomcellstates[greenviolinplotandthedashedlineinFig.

2(d)

].Bysettingthestateofeachgenetobe0or1randomly,wegeneratedmultiplerandomcellsstatesandobtainedtheirfrustrationscoredistribution.Thefrustrationscoreofthecellularsystemislowerintheinitialand?nalmetastablestatescomparedtothatoftherandomstates,butitreachesapeakvaluethatishigherduringthetransition.Theseresultsalignwithourexpectationthatthetransitionsarecoordinatedratherthanrandom.

Duringtransition,cellsneedtoexploreregionsofthestatespacerarelyexploredbycellsinastablephenotype,lead-ingtoincreasedcellstateheterogeneityduringtheprocess

[7]

.Weusedadistributionofstateoverlaptoquantifytheheterogeneityduringcellphenotypictransition.OriginallyusedinspinglassandBooleanstatemodelstocharacterizestatestructures[

8

,41],thestateoverlapisde?nedas

qab=

sampledinthesameVoronoigridcorrespondingtoareactioncoordinatepointi.Nisthenumberofgenes,andthesumisoverallthegenesincludedintheanalyses.Notethatqab=1

ifthestatevectorss=s,andqab=?1ifss.Lower

heterogeneitygivesanarrowerdistributionandahighermeanvalueofthestateoverlap.Figure

2(e)

showsthedistributionofstateoverlapPr(qab)calculatedateachreactioncoordinatepoint.Duringtransition,thedistributionalongthereactioncoordinateinitiallybecomesmoredispersive,anditsmeanvaluegraduallyapproaches0from1.Then,thedispersionnarrowsdownandreturnstoameanvaluecloseto1.Thistrendisapparentfromthevariationofthefullwidthathalf-maximum(FWHM)ofthePr(qab)alongthereactionco-ordinate[Fig.

2(e)

inserted].Thisre?ectsatemporalincreaseinheterogeneityduringthetransition.

Wealsoquanti?edgeneregulatorynetworkheterogeneity

[42,43].Networkheterogeneitymeasureshowhomogenously

theconnectionsaredistributedamongthegenes[

42]

.Highheterogeneityindicatesthecoexistenceofhubgeneswithhighconnectivityandgeneswithlowconnectivity.Thenet-workheterogeneityinitiallyincreases,reachingamaximuminneuroblastcells,andthendecreasesasthecellsapproachthematuregranulecellstate(Fig.S1b)[

42]

.Weobservedasimilarpatternusinganothertypeofheterogeneitycalleddegreeheterogeneity,whereahighervaluere?ectsamoreunevendistributionofthedegreeofgenesinageneregulatorynetwork.Thedegreeofageneisde?nedasthenumberof

connectionsthisgenehastoothergenes(Fig.S1c)[

43]

.Thisresultfurtherindicatesthattheobservedtrendisnotduetoaspeci?cchoiceofheterogeneitymeasure.

D.Reconstructedgeneregulatorynetworkreveals

increasedintercommunityinteractionsatanintermediate

stageoftransition

Toexaminethenatureoftheincreasedinteractions,wedividedtheinferredgeneregulatorynetworkintofourcom-munitiesusingtheLeidenmethod(MaterialsandMethodsintheSM[

29

])[

44

,45]

.Thenumberofeffectiveintracom-munityedgescorrelateswiththenumberofactivegenes(Fig.S1d).Wedidnotobserveauniversalpatternamongthefourcommunitiesregardinghowintracommunityinteractionschangealongthereactioncoordinate.Speci?cally,thenum-berofintracommunityedgesforcommunity0increases,forcommunity2decreases,whilethoseforcommunities1and3peakinthemiddle(Fig.S1dleft).Incontrast,thenumberofeffectiveintercommunityedges?rstincreasesandthendecreasesbetweenpairsofallfourcommunities[Fig.

2(f)]

.Theintercommunityinteractionstrengths,de?nedasthetotalnumberofeffectiveedgesbetweendifferentcommunities,arethestrongestatr=5[Fig.

2(g)]

.Thisvariationininter-communityeffectiveedgesdoesnotcorrelatewiththetotalnumberofactivegenes,withacorrelationcoef?cientonly0.26(Fig.S1e).

E.Theobservedpatternsarerobustwithdifferent

transitionmodels

Notethatareactioncoordinatere?ectsthesequenceofeventsduringthetransitions,primarilyusingthedirectionalinformationofcellstatetransitionswithoutreferringtotheactualratesofeachtransitionstep.Consequently,weexpectthatourconclusionshereareinsensitivetotheassumptionβi=1whileestimatingthesplicing-basedRNAvelocities.Tosupportthisexpectation,werepeatedtheaboveanalysesusinganotherdiscretetransitionmodelfromRNAvelocitiesbasedonacosinecorrelationkernelmethodofBergenetal.

[23],

whichiswidelyusedintheRNAvelocity?eld.Mathematicalanalyseshaveshownthatthecosinecorrelationkernelmethodasymptoticallykeepsthedirectionalinformationofaveloc-ityvectorbutcompletelylosesthemagnitudeinformation

[46]

.Theanalysesindeedproducedareactioncoordinate,andthepatternsofthefrustrationscoreandthenumberofeffectiveintercommunityedgeschangedalongthereactioncoordinatesimilartowhatwasobservedwithGraphDynamo(Fig.S2a-c).

Weperformedtheabovenetworkanalysesusing678bi-narizedgenes.GiventhatanscRNA-seqdatasetcantypicallybewellrepresentedasalow-dimensional(e.g.,<30)manifoldembeddedinthehigher-dimensionalgenespace,weexpectedthatthegenesselectedthroughtheaboveproceduresfaithfullyrepresentthedatamanifoldbasedonthecelebratedWhitneyembeddingtheorem[

47]

.Tofurtherruleoutthein?uenceofgeneselectioninoursubsequentanalyses,wealsoadoptedanadditionalproceduretodeterminethethresholdforbinarizinggenesthatdonotshowswitchlikebehavior(MaterialsandMethodsintheSM[

29

]).Forthedentategyrusneurogenesis

WANG,NI,POE,ANDXINGPRXLIFE2,043009(2024)

043009-6

FIG.3.Analysesonhematopoiesisdatasetrevealsfeaturesofaconcertedmechanism.(a)TransitiongraphofhematopoiesisbasedonRNAvelocity.Colorrepresentscelltype.Arrowindicatesthedirectionoftransition.(b)Reactioncoordinate(largecoloreddots,startsfromblueandendsinred)ofhematopoiesiswithcorrespondingVoronoigrids.Smalldotsrepresentcellswithcolorsindicatingcelltypes.

(c)Frustrationscorealongthereactioncoordinateinhematopoiesis.(d)Variationofdistributionofstateoverlapalongthereactioncoordinateinhematopoiesis.Colorsofthedistributionscorrespondtothatofreactioncoordinatepoints.Thedistributionsofinitialand?nalreactioncoordinateareplottedwithadottedline.InsertedgraphistheFWHMofthedistributionsalongthereactioncoordinate.(e)Cell-speci?cvariationofeffectiveintercommunityregulationinhematopoiesis.Colorrepresentsthenumberofeffectiveintercommunityedgeswithineachcellinthegeneregulatorynetwork.(f)Evolutionofthenumberofeffectiveintercommunityedgesalongthereactioncoordinateduringhematopoiesis.Eachnoderepresentsagene.Thecolorofthenodeindicatestheindexofthecommunity.Arrowsrepresentdirectionofregulation.ristheindexinthereactioncoordinate.

dataset,thisprocedurebinarizedallthe1960high-variancegenes.Theresultscomputedwiththesegenes(Fig.S3a)aresimilartothoseinFig.

2.

Insummary,theabovenetworkanalysessupportthecon-certedmechanism[Fig.

2(h)]

.Coexpressionofcon?ictinggenesleadstoincreasedintercommunityedges,andfrustratededges.Someofthegenestransientlyactashubgenes,leadingtoincreasednetworkheterogeneity.

F.R