StatisticalPhysics of the SymmetricGroup Mobolaji Williams Department of Physics, Harvard University, Cambridge, MA 02138, USA ∗ Eugene Shakhnovich DepartmentofChemistryandChemicalBiology,HarvardUniversity,Cambridge,MA02138, USA† (Dated: March7,2017) 7 1 Orderedchains(suchaschainsofaminoacids)areubiquitousinbiologicalcells,andthesechains 0 performspecificfunctionscontingentonthesequenceoftheircomponents. Usingtheexistenceand 2 generalpropertiesofsuchsequencesasatheoreticalmotivation,westudythestatisticalphysicsofsys- temswhosestatespaceisdefinedbythepossiblepermutationsofanorderedlist,i.e.,thesymmetric r a group,andwhoseenergyisafunctionofhowcertainpermutationsdeviatefromsomechosencorrect M ordering.Suchanon-factorizablestatespaceisquitedifferentfromthestatespacestypicallyconsid- eredinstatisticalphysicssystemsand consequentlyhas novelbehavior insystemswithinteracting 5 andevennon-interactingHamiltonians. Variousparameterchoicesofamean-fieldmodelrevealthe systemtocontainfivedifferentphysicalregimesdefinedbytwotransitiontemperatures,atriplepoint, ] andaquadruplepoint.Finally,weconcludebydiscussinghowthegeneralanalysiscanbeextendedto h statespaceswithmorecomplexcombinatorialpropertiesandtootherstandardquestionsofstatistical c mechanicsmodels. e m - I. INTRODUCTION t a t s Chains of amino acids are important components of . t biological cells, and for such chains the specific order- a ingoftheaminoacidsisoftensofundamentaltothere- m sulting functionandstabilityof the foldedchainthatif - majordeviationsfromthecorrectorderingweretooccur, d n thefinalchaincouldfailtoperformitsrequisitefunction o withinthecell,provingfataltotheorganism. FIG.1: Foldingvs. Design(orInverseFolding) c More specifically, we see the relevance of correct or- problems: Theproteinfoldingproblemisconcerned [ dering in the study of protein structure, which is of- withdeterminingthethreedimensionalstructure 2 tendividedintothe proteinfoldingandproteindesign producedbyaparticularsequenceofaminoacids. The v problem. While the protein folding problem concerns proteindesignproblem(whichmotivatesthecurrent 4 findingthethree-dimensionalstructureassociatedwith work)isconcernedwithfindingthesequence(s)of 7 agivenaminoacidsequence,theproteindesignproblem aminoacidswhichyieldagiventhreedimensional 3 (alsotermed the inverse-foldingproblem; see Figure 1) polypeptidestructure. Anumberofapproachestothe 7 concerns finding the correctamino acidsequence asso- designproblemaregivenin[1] 0 ciatedwithagivenproteinstructure. . 1 Anaspectofonesolutiontotheproteindesignprob- 0 7 lem is to maximize the energy difference between the For example, for a polypeptide chain with N residues, 1 low-energy folded native structure and the higher en- ratherthansearchingovertheentiresequencespace(of : ergymisfolded/denaturedstructures. Indoing so, one size20N),onesearchesoveraspaceofsequences(ofsize v takes native structure asfixed and then determines the N!/n !n !...n !)whicharedefinedbyafixednumber i 1 2 20 X sequence yielding the minimum energy, under the as- ofeachaminoacid. sumption (termed the ”fixed amino-acid composition” r This aspect of the protein design problem alerts one a assumption) thatonlycertainquantitiesof amino-acids to agap in the statisticalmechanics literature. Namely, appear in the chain [2]. In this resolution (specifically theredonotseemtobeanysimpleandanalyticallysolu- termedheteropolymermodels[3][4])thecorrectamino blestatisticalmechanicsmodelswherethespaceofstates acid sequence is found by implementing an MC algo- isdefinedbypermutationsofalistofcomponents. rithminsequencespacegivenacertainfixedaminoacid We can take steps toward constructing such a model composition. This entailsassuming the number of var- by considering reasonable general properties it should ioustypesofaminoacidsdoesnotchange,anddistinct have. Ifweassumetherewasaspecificsequenceofcom- statesinsequencespacearepermutationsofoneanother. ponentswhichdefinedthelowestenergysequenceand wasthermodynamicallystableinthemodel,thendevi- ationsfromthissequencewouldbelessstable. Because ∗[email protected] oftherolesequencesofmoleculesplayinbiologicalsys- †[email protected] tems, it is worth asking what features we expect such 2 sequences to have from the perspective of modeling in pressedincomponentformas statisticalmechanics. InSectionIIweintroducethemodel,andcomputean ω~ =(ω1,...,ωN). (2) exact partition function which displays what we term Wewilltakethecollectionofcomponents ω asintrin- “quasi”-phasetransitions—atransitioninwhichthese- { k} sictooursystem,andthustakethestatespaceofoursys- quence of lowest energy becomes entropically disfa- temtobethesetofallthevectorswhoseorderingofcom- vored above a certain temperature. In Section III, we ponentscanbeobtainedbypermutingthecomponents extendthepreviousmodelbyaddingaquadraticmean of ~ω, i.e., all permutations of ω ,...,ω . We represent field interaction term and show that the resulting sys- 1 N temdisplaystwotransitiontemperatures,atriplepoint, anarbitrarystateinthisstatespaceas~θ = (θ1,...,θN), andaquadruplepoint. InSectionIV,wediscussvarious wheretheθk aredrawnwithoutrepeatfrom ωk . For- { } ways we can extend this model in theoretical or more mally,wewouldsayourspaceofstatesisisomorphicto phenomenologicaldirections. thesymmetricgroupon~ω([5]). Wewillthusdenoteour statespaceas Sym(ω) := SetofAllPermutationsof(ω ,...,ω ). 1 N II. SYSTEMANDPARTITIONFUNCTION (3) andthen anarbitrarystate~θ isjustanelementelement Ourlargergoalistostudyequilibriumthermodynam- ofthisset. ics for a system defined by permutations of a set of N Asafirstformulationofthemodel,wewilltake~θ0 =ω~ componentswhereeachunique permutationisdefined (thecorrectsequence)torepresentthezeroenergystate byaspecificenergy. Ingeneral,weshouldconsiderthe inthesystem,andforeachcomponentθiofanarbitrary casewherethesetofN componentsconsistsofLtypes vector ~θ which differs from the corresponding compo- ofcomponentsforwhichifnkisthenumberofrepeated nentωiin~ω,thereisanenergycostofλi >0. TheHamil- componentsoftypek,then L n = N. Forsimplic- tonianisthen k=1 k ity, however, we will take n = 1 for all k so that each kP N componentisofauniquetypeandL=N. ( θ )= λ I , (4) To study the equilibrium thermodynamics of such a HN { i} i θi6=ωi i=1 systemwithafixedN atafixedtemperatureT,weneed X to compute its partition function. For example, for a where θ and ω are components of vectors ~θ and ω~ re- i i sequence with N components (withno components re- spectively,andI isdefinedby peated),thereareN!microstatesthesystemcanoccupy andassumingwelabeleachstatek = 1,...,N! 1,N!, 1 ifAistrue andassociateanenergyǫkwitheachstate,thent−hepar- IA ≡ 0 ifAisfalse. (5) titionfunctionwouldbe (cid:26) N! Wenotethatalthoughwelabelourgeneralstateasθ~ = Z = e−βǫk, (1) (θ ,...,θ ),thecomponentsθ ,...,θ canonlytakeon 1 N 1 N kX=1 mutually-exclusivevaluesfromtheset{ωk}. Wewanttoexploretheequilibriumthermodynamics whereǫ foreachstatekcouldbereasonedfromamore k ofasystemwithaHamiltonianofEq.(4). Thisamounts precisemicroscopictheoryofhowthecomponentsinter- tocalculatingthepartitionfunction actwithoneanother. Phenomenologically,Eq.(1)would be the most precise way to construct a model to study N theequilibriumpropertiesofpermutations,butbecause Z ( βλ )= exp β λ I , (6) itbearsnoclearmathematicalstructure,itisunenlight- N { i} ~θ∈SXym(ω) − Xi=1 i θi6=ωi! eningfromatheoreticalperspective. Instead we will postulate a less precise, but theoreti- whereSym(ω)isagainthesetofallpermutationsofthe cally more interesting model. For most ordered chains components of (ω ,...,ω ). To find a closed form ex- 1 N in biological cells, there is a single sequence of com- pression for the partition function, we group terms in ponents which is the “correct” sequence for a partic- Eq.(6) according to the number of ways to completely ular macrostructure. Deviations from this correct se- reorderjcomponentsin~ωwhilekeepingtheremaining quenceareoftendisfavoredbecausetheyformlessstable componentsfixed. Eachsuchreordering(i.e.,eachvalue macrostructuresortheyfailtoperformtheoriginalfunc- of j) is associated with a sum over products of e−βλi tionofthe“correct”sequence. Withthegeneralproper- termswithjfactorsofe−βλi (forvariousi)ineachterm. tiesofsuchsequencesinmind,wewillabstractlyrepre- Thetotalpartitionfunctionisasumofallsuchreorder- sentoursystemasconsistingofN siteswhicharefilled ingsforalljsfrom0toN. Ascanbeseenfromadirect with particular coordinate values denoted by ω . That k is,wehaveanarbitrarybutfixedcoordinatevector~ωex- 3 expansionofEq.(6),wehave definedas N N g (j)= d (13) Z ( βλ )= d Π e−βλ1,...,e−βλN , (7) N j j N { i} j j (cid:18) (cid:19) j=0 X (cid:0) (cid:1) isthenumberofwaystoreorderalistofN elementsso whered ,termedthenumberofderangementsofalistof that j elements are no longer in their original position. j j([6]),isthenumberofwaystocompletelyreorderalist Thiscombinatorialdefinitionofg (j)willproveuseful N of j elements. The quantity Π (x ,...,x ), termedthe when we explore the phase behavior of more complex j 1 N jth elementary symmetric polynomial on n ([7]), is the modelsofpermutations. sumofallwaystomultiplyjelementsoutoftheN term FromtheformofEq.(10),itisclearthat,physically,its set x1,...,xN . For example, Π2(x1,x2,x3) = x1x2 + associated Hamiltonian is not realistic as it places dis- { } x2x3+x3x1. Bydefinition tinct permutations (which in any true physical system mostlikelyhavequitedifferentenergyproperties)inthe 1 dk N same degenerate energy state. Still, from a theoretical Π (x ,...,x )= (1+qx ) . (8) k 1 N k! dqk i perspective,thesimplicityofthismodelmakesitasuit- " # Yi=1 q=0 ablestartingpointforstudyingthegeneralpropertiesof systemsofpermutations. By the definition of the incomplete gamma function as Γ(x,α)= ∞dttx−1e−tanditsrelationtoderangements α (i.e.,d =Γ(j+1, 1)/e,see[8]),wethenfind j R − A. Phase-likeBehaviorof“Non-Interacting”System ∞ N Z ( βλ )=e−1 dte−t tjΠ e−βλ1,...,e−βλN N i j { } Wecaninvestigatethephase-likebehaviorofthesys- Z−1 j=0 X (cid:0) (cid:1) tem defined by the Hamiltonian Eq.(4) (for constant λi ∞ N across i), byapplying steepestdescentto Eq.(11) in the = dse−s 1+(s 1)e−βλℓ , (9) − N 1limit. Doingso,wefindwecanapproximatethe Z0 ℓY=1h i fre≫eenergyofthesystemtobe whichisthedesiredclosed-formexpressionforthepar- βF = lnZ (βλ ) titionfunctionofthissystem. − N 0 WithEq.(9),theproblemofabstractlystudyingather- Nβλ eβλ0 N 1 +F (N), (14) 0 0 ≃ − − − malsystemofpermutationswithHamiltonianEq.(4)is, from the perspective of equilibrium statistical mechan- whereF0(N)isaβλ0in(cid:0)dependentcon(cid:1)stant. Notingthat ics,nowcomplete. However,therearestillsomephysi- j = ∂lnZN(βλ0)/∂(βλ0) = ∂F/∂λ0,wefindtheav- h i − calandtheoreticalresultswhichcanbeteasedfromthis erage number of incorrect components satisfies the fol- formalism. Specifically,wecanaskwhetherthissystem lowingequationofstate: exhibits phase transitions. To answer this question, it wouldprovemoreanalyticallytractabletotakeλ = λ j N eβλ0. (15) i 0 h i≃ − foralli. Withthiscondition,ourpartitionfunctionsim- ByEq.(12),wecaninferthat j mustbegreaterthanor plifiesto h i equal to 0. However, the right-hand-side of Eq.(15) ex- N hibitsnosuchexplicitconstraint. Thuswecaninferthere Z (βλ )= g (j)e−jβλ0 (10) is a phase-like transition in our system at the tempera- N 0 N ture j=0 X ∞ = dse−s 1+(s 1)e−βλ0 N (11) k T = λ0 . (16) Z0 − B c lnN (cid:0) (cid:1) wherewetransformedourHamiltonianas N( θi ) Belowthistemperature,wemusthave j 0andthus H { } → h i ≃ (j)=λ0jwithj,definedas the“correctpermutation”hasthelowestfreeenergyand H is thermodynamically favored; above this temperature, N j > 0 and the system is in a disordered phase where j I , (12) h i ≡ θi6=ωi thepreviouslowestenergy“correct-permutation”isen- Xi=1 ergeticallydisfavored. Interestingly, this transition arises from the naively the number of components of ~θ which are not equal to non-interactingHamiltonian the correspondingcomponent in~ω. Wecallj the num- ber of incorrect components of ~θ, and if j = N we say N ~θ is completely disordered. The factor gN(j) in Eq.(10) is HN({θi})=λ0 Iθi6=ωi. (17) i=1 X 4 1. NotaTruePhaseTransition Weclaimthesystemdoesnotexhibittruephasetran- sition behavior because many of these results are not consistent with the traditional thermodynamic defini- tionofphasetransitions. For one, phasetransitionsare associatedwithdivergencesinthederivativesofthefree energy, but there is no divergence in the free energy associated with the partition function Eq.(11) for pos- sible parameter values. Also, the result Eq.(15) arises from the steepest descent approximation which makes j ’stemperaturedependencenear j = 0appearnon- h i h i FIG.2: Freeenergyfor“Non-InteractingModel”: For differentiablewhenbyEq.(11)itisactuallydifferentiable λ0 >0,theLandaufreeenergyofthesystemasa overitsentire domain. Finally, withEq.(10) wecande- functionofj,thenumberofincorrectcomponentsin~θ, fineaLandaufreeenergyF(j)forthissystemaccording isalwaysconvexwithasingleglobalminimum. toZ = e−βF(j),andwhatwemayordinarilylabelas j Becausej 0,thej <0domainofeachplot(dashed a phase transition (i.e., going from j = 0 to j = 0) ≥ P h i h i 6 section)isinaccessible. Forsufficientlylow arises, not from changes in the functional form of the temperatures,theminimumisatj =0,butaswe Landaufree energy as we see in real phase transitions, increasethetemperaturebeyondEq.(16),thefree butfromchangesintheexcludedregionoftheLandau energycurvemovestotheright(butretainsits freeenergy(SeeFigure 2). Becausethe functionalform functionalform)andthenewminimumisataj >0 ofthefreeenergyremainsthesameweobservenotrue value. ActualplotsofEq.(20)forj <0requireusto phasetransition. replacethecombinatorialargumentofthelogarithm Alternatively, a heuristic argument for the non- withitscorrespondinggammafunctionexpression. existenceofphasetransitionsinourpermutationmodel ismathematicallyverysimilar totheLandauargument ([9])forthenon-existenceoftransitionsin1dIsingMod- els. For our permutation system with N lattice sites, thestateofzeroenergyandzeroentropyconsistsofev- erysitebeingoccupiedbyitscorrectcomponent. Toin- creasetheenergyofthissystem,wecanchoosejsitesto containincorrectcomponents,thusgivingusanenergy We say “naively non-interacting” because Eq.(17) con- = λ j. The number of wayswe can choose these j j 0 H sists of a sum over linear functions of a single index i, componentsisgivenbyEq.(13)Thus,uponintroducing andthus doesnot suggestanycoupling betweenterms j = 0 incorrect components, the change in the Landau 6 ofdifferingindex. However,statisticalmechanicstellsus freeenergyofoursystemis thattheenergyofasystemisn’ttheonlythingwhichde- termines the thermodynamic behavior of a system. In- N ∆F(j)=λ j k T ln d j(λ k T lnN), deedwehavetoconsiderentropiccontributionsaswell, 0 − B j j ≃ 0− B (cid:20)(cid:18) (cid:19) (cid:21) and in this system the entropy (as it is a function of j) (20) candrivethermodynamic behavior. Inotherwords, al- though the Hamiltonian is depicted as non-interacting where we took these results in the 1 j N limit ≪ ≪ andcanset-theoreticallyberepresentedas and used d j!/e. In the thermodynamic (N ) j ≃ → ∞ limit, we find that ∆F(j) meaning there is no system = 1 2 N, (18) non-zeroT atwhichtheen→trop−i∞ccontributionbecomes H H ⊕H ⊕···⊕H subdominanttotheenergy. Thusthesystemexhibitsno our systemreallyexhibits interactions betweencompo- phasetransition. nents because our total space of states cannot be fac- S torized: = ... . (19) system 1 2 N S 6 S ⊗S ⊗ ⊗S III. PARTITIONFUNCTIONFORINTERACTING Thus the “non-interacting” system exhibits a transition MODEL atEq.(16)duetothecouplednatureofthestatespace. As discussedinthesubsequentsection,wetermthistransi- tiona“quasi-phasetransition”becauseitdoesnotbear Whenwefirstconsideredamodelofthermalpermu- allofthestandardpropertiesweexpectinphasetransi- tations,webeganwithaHamiltonianwheresitesdidnot tions. interactwithoneanotherandeachhadasite-dependent 5 energycostforbeingincorrectlyoccupied: A. CalculatingOrderParameter ( θ )= λ I . (21) H { i} i θi6=ωi i X More generally, we can consider Hamiltonians with an arbitrary number of multiple-site interaction terms. Our goal is to analyze the “quasi”-phase behavior of SuchaHamiltoniancouldbewrittenas this system in a way analogous to our analysis for the non-interactingsystem. TodosowebeginwiththeLan- 1 ( θ )= λ I + µ I I + . daufreeenergyfunction H { i} i θi6=ωi 2 ij θi6=ωi θi6=ωi ··· i i,j X X λ 1 (22) F (j,β)=λ j+ 2 j2 lng (j). (26) N 1 N The first term in Eq.(22) associates an energy cost of λ 2N − β i with incorrectly occupying the component at position Alternative starting points for this derivation are pre- i. The second term models interactions between sites sentedinAppendixA.Oursystemisconstitutivelydis- where the correct (or incorrect) occupation of a single crete,soitisnotpreciselycorrecttodiscussourfreeen- site determines the energy of another. The exact val- ergy in the language of analysis, but given our expres- ues of λ and µ could be chosen to ensure the “cor- i ij sionforEq.(13)wecanmapthissystemtoacontinuous rect”state(θ = ω foralli)isnon-degenerateasinthe i i one which bears the same thermodynamic properties non-interactingmodel. Theellipsesrepresenthigheror- andforwhichanalysisisappropriate. Specifically,ifwe derinteractionsinthisframework.Hamiltonianssuchas takejtobecontinuousandusetheidentityΓ(x+1)=x!, Eq.(22)shouldbemorephysicallyrelevantastheywould wecanwrite correspond to systems where the energy cost for devi- atingfromthelowestenergypermutationisnotsimply Γ(N +1) Γ(j+1, 1) linearbutcouldberepresentedasatensorvaluedfitting gN(j)= − , (27) Γ(j+1)Γ(N j+1) e function. − WiththeapproximationΓ(j+1, 1) Γ(j+1)andthe − ≃ substitutionEq.(27),Eq.(26)thenbecomes λ 1 Wecanmakeprogressinstudyingthethermodynam- f (j,β)=λ j+ 2 j2+ lnΓ(N j+1)+f (28) N 1 0 icsofmoregeneralHamiltonianslikeEq.(22)byfirstonly 2N β − considering first- and second-order interaction terms where we defined our approximated free energy as andtakingthe interactionstobeconstants: λ = λ for i 1 f (j,β) and collectedthe j independentconstants into alli;µ = λ /N foralli,j. Thefactorof1/N ischosen N ij 2 f . Now Eq.(28) is fully continuous and amenable to sothatthesecondtermmatchestheextensivescalingof 0 analysis. To find the thermodynamic equilibrium of thefirstterm. Thepartitionfunctionforsuchparameter this system, we need to find the value of j for which selectionsisthen ∂f (j,β)/∂j = 0 and ∂2f (j,β)/∂j2 > 0. For the first N N N conditionwehave Z (β;λ ,λ )= exp βλ I N 1 2 − 1 θi6=ωi− ∂ λ2 1 ~θ∈SXym(ω) Xi=1 ∂jfN(j,β)=λ1+ Nj− βψ0(N −j+1)=0. (29) N βλ 2 I I , (23) Forx 0.6wehave 2N θi6=ωi θj6=ωj ≥ i,j=1 X ψ (x) ln(x 1/2), (30)  0 ≃ − whereλ andλ areinteractionparameterswithunitsof 1 2 ascanbeaffirmed byTaylor expansionor plotsof each energy. Wecanalsowrite this partition functionin the side. Since the argument of ψ (N j +1) is bounded Eq.(12)basisas 0 − belowby1,theapproximationinEq.(30)canbeapplied N toEq.(29). Withthissubstitution, andsettingtheresult Z (β;λ ,λ )= g (j)e−βE(j), (24) tobevalidfortheequilibriumvaluej =j,wethenfind N 1 2 N theconstraint j=0 X whereg (j)isdefinedinEq.(13)and eβλ2j/N = e−βλ1 j N 1/2 , (31) N − − − λ whichhasthesolution (cid:0) (cid:1) (j)=λ j+ 2 j2 (25) 1 E 2N j 1 βλ 1 =1 W 2eβλ1+βλ2 + , (32) istheenergyfunctionforthesystem. N − βλ N O N 2 (cid:18) (cid:19) (cid:18) (cid:19) 6 where W is the (branch unspecified) Lambert W func- tion,definedby W(xex)=x. (33) Tospecifythe branchof the W which corresponds to a stableequilibriumwecomputethesecondderivativeof our free energy at this derived critical point. Doing so yields ∂2 1 1 1 f (j =j,β) λ + ∂j2 N ≃ N 2 β1 j/N (cid:18) − (cid:19) λ 1 2 = 1+ . (34) N  W βλ2eβλ1+βλ2  N  (cid:16) (cid:17) ThusEq.(32)(forλ > 0)yieldsafreeenergyminimum 2 for FIG.3: (Coloronline)Possiblefunctionalformsof Eq.(28): Wenotethatthestabilitiesthatdefinethej =0 βλ W 2eβλ1+βλ2 > 1, (35) andj =N pointsarenotthermodynamicstabilities (cid:18) N (cid:19) − (namelytheydon’tarisefromthef′(j)=0condition). Rathersincethespectrumofj valuesisboundedbelow and yields a maximum for the inverse condition. This by0andabovebyN,owingtotheseboundary amountstostatingthatthestableequilibriumforjisde- conditionsthesystemcanbecometrappedinordinarily finedbytheprincipalbranchoftheLambertWfunction unstablepartsofthefreeenergycurve. Thecolors whereW =W 1,andtheunstableequilibriumforj 0 ≥− matchthecoloroftheassociatedregionofparameter isdefinedbythenegativebranchwhereW =W < 1. −1 − spaceinFigure4. Thus,theorderparameterforthissystemis j 1 βλ 1 0 =1 W 2eβλ1+βλ2 + . (36) tionsinthissystem,thereareinfactfivedistinctregimes 0 N − βλ2 (cid:18) N (cid:19) O(cid:18)N(cid:19) ofparameterspace. Wecanobtainaqualitativesenseof theseregimesbycreatingschematicplotsofthefreeen- Wenotethattakingλ 0andusingW(x)=x+ (x2) 2 → O ergyEq.(28)for variousparametervaluesof λ1 and λ2. for x 1 returns us to the non-interacting result | | ≪ Thepossibleplotscanbeplacedintofivecategoriesac- Eq.(15). cording to the plot’s stable or metastable j values. We For completeness, we define the value of j which depictthesepossibleplotsinFigure3. Wenotethatonly yields a free energy maximum as j ; it is related to −1 thefreeenergyplotswithvalidvaluesofj containwhat 0 Eq.(36) by replacing the principal branch function W 0 wenormallyconsiderathermodynamicequilibrium;the withW . −1 otherplotshave“stable”valuesofjarisingonlyfromthe j =0and/orj =N boundaryconditions. Qualitatively,wecannamethestatesaccordingtothe B. DiscussionofParameterSpace sequence space to which their equilibrium values of j correspond. We know for j = 0, our system is in a In the previous section, we found that the order pa- statewithzeroincorrectcomponentsinθ~andhencethe rameterforthissystemwasgivenbyEq.(36). Wenoted system is “perfectly ordered” or just ”ordered”. Con- that this solution represented a local minimum of the versely for j = N our system has N incorrect compo- freeenergyaslongastheLambertWfunctionsatisfied nentsandhence the systemis“completelydisordered” W = W0 > 1. Thus when this condition is violated, orjust“disordered”.Thein-betweencaseofj =jwhere − j is no longer a valid stable equilibrium, and our sys- 0 < j < N can be given the relatedlabel of “partially- 0 temhasundergonea“quasi”-phasetransitionorsimply ordered”. Thus, the regime names associated with our atransition. possiblevaluesoftheorderparametersare Moreover,our valuesofj areboundedbelowbyj = OrderedRegime(j = 0stable):Neitherj orj 0 and bounded above by j = N, neither conditions of • 0 −1 exist;f(N,β)>0. which are naturallyconstrained byEq.(36). Thus these twoconditionsareassociatedwithtwoothertransitions. DisorderedRegime(j = N stable):Neitherj or Inall,then, therearethreeconditions whichdefinethe • j exist;f(N,β)<0. 0 −1 quasi-phaseboundariesinthissystem. While therearethreeconditions which definetransi- Partially-Ordered Regime (j = j0 stable): Only • 7 TABLEI:Functionsdefiningboundariesbetween ischaracterizedbythecondition parameterregimes. λ =lnN/β and λ = 1/β. (37) 1 2 − RegimeTransition BoundaryinParameterSpace Theseconditionscharacterizethesystem’striplepoint. 1 Order/Partial-Order λ1 = lnN Similarly, for N 1, the Partially-Ordered, Dis- β ≫ ordered, Order-Partially-Ordered Metastability, and 1 N Order-Disorder Metastability coexistence point is char- Order/Order-Partial- λ = W e−βλ1 2 β −1 − e acterizedbythecondition OrderMetastability (cid:18) (cid:19) λ Order/Disorder λ2 =−1 11/N λ1 =lnN/β ≃−λ2. (38) − Thisconditioncharacterizesthe quadruplepointof the system. j exists. Figure4alsodepictsthe possible regimesof oursys- 0 tem for a given temperature and various Hamiltonian Order and Disorder Metastable Regime (j = 0 parametersλ1andλ2. Morephysically,wemaybeinter- • andj = N stable): Onlyj exists. estedinknowingwhatarethe“quasi”phaseproperties −1 ofasystemwithafixedλ andλ andavariabletemper- 1 2 Order and Partial-Order Metastable Regime ature. That is, what are the temperatureswhich define • (j = 0andj = j0stable): Bothj0andj−1exist. thevarioustransitionsbetweenregimesinthesystem? AnarbitrarypermutationsystemforafixedN andat Wenote thatitseemstobe afundamentalfeature(or a avariabletemperatureischaracterizedbyaspecificen- lack of one) of this system, that the free energy Eq.(28) ergyfunctionEq.(25). Suchasystemisthereforedefined doesnotadmitametastabilitybetweenpartial-orderand byaspecificλ andλ ,andthesystemcanbeassociated 1 2 disorder. withaparticularpoint(andhenceregion)intheparame- terspaceofFigure4. Aswevarythetemperatureofthis system, the temperature dependent regime-coexistence C. Monte-CarloGeneratedParameterSpace lines change and if theychange in such awayasto ex- tendtheregionofaregimetonewlyencompassourorig- inalpointthenoursystemhasundergoneatransition. In With these regime definitions, we can depict the pa- thisway,wecandefinethetemperatureswhichcharac- rameter space graphically. In Figure 3 we showed the terizevariouspossibletransitionsofthissystem. possible formsofthe freeenergyforthissystemwhere First,fromFigure4andEq.(C10)wenotetheregime- eachwascategorizedaccordingtotheexistenceofthelo- coexistence line between the partially-ordered and dis- calminimum criticalpointj , the existenceof the local 0 ordered regime is independent of temperature, and so maximumcriticalpointj ,andthesignofthequantity −1 there is no criticaltemperature defining a partial-order f (j,β). Wecanextrapolatethiscategorizationtoλ λ N 1 2 − todisordertransition. parameter space, by determining which regions of pa- FromEq.(C3),wecaninferthatthepartial-ordertoor- rameterspacecorrespondtospecificplotsinEq.(28). Do- dertransitionischaracterizedbymovingbelowthetem- ingsothroughtheMonteCarloproceduredescribedin perature Appendix B, we generated 10,000points of the param- eter space diagram in Figure 4 for β set to 1. We note λ thattheparameterspaceexhibitsfiveregimesseparated k T = 1 . (39) B c1 lnN bythreelinescitedinTableIeachofwhichcorrespond to the three conditions mentioned at the beginning of Contingentonwhichregionofparameterspacethesys- this section. These lines can be derived analytically(as tem lies, this temperature also characterizes the disor- showninAppendixC)byconsideringtheconditionsin der to order-disorder metastability transition and the turnandwhichregimestheyservetoconnect. partial-ordertoorder-partial-ordermetastabilitytransi- tion. AndfromEq.(C6),wecansolvefortheassociatedtran- D. TripleandQuadruplePointsandTransition sitiontemperaturegivenfixedλ andλ tofind 1 2 Temperatures N −1 k T =(λ +λ ) W (λ +λ ) (40) From Figure 4 we see that our system is character- B c2 1 2 0 −eλ 1 2 (cid:20) (cid:18) 2 (cid:19)(cid:21) ized by two points where there is a coexistence be- tween multiple regimes. Given that W ( e−1) = wherethisexpressionisonlyrelevantfor λ < λ < 0 −1 1 2 − − 1, we have that the Ordered, Partially-Ordered, and and λ > k T lnN. Moving above this temperature 1 B − Order-Partially-OrderedMetastabilitycoexistencepoint leads to the order to order-partial-order metastability 8 FIG.4: (Coloronline)λ λ ParameterSpaceforInteractingMeanFieldSystem: Wesetβ =1andN =100. We 1 2 − followedtheMonteCarloprocedureoutlinedinthenotesfor10,000points. Inthefigurewealsodenotedthe analyticlines(i.e.,Eq.(C3),Eq.(C7),andEq.(C10))whichdefinetheseparationbetweenthephases. Thecolors correspondtothecolorsofthefreeenergycurveinFigure3. transition. of typical canonical models in statistical mechanics. In particularwehopetoextendittonon-trivialsitedepen- dentinteractions. Forexample,anearestneighborinter- action Hamiltonian of the kind which characterize the IV. DISCUSSION IsingModel, In this work, motivated by an abstraction of a foun- N dationalprobleminproteindesign,wepositedandana- ( θ )= q I I , (41) lyzedthebasicpropertiesofastatisticalphysicsmodelof H { i} − θi6=ωi θi+16=ωi+1 i=1 X permutations. Formally,weconsideredasimplestatisti- calphysicsmodelwherethespaceofstatesforN lattice would be an alternative physical extreme to the mean- siteswasisomorphic tothe symmetric group ofdegree fieldinteractionsconsideredinSectionIII. N [5], and where the energy of each permutation was We could also consider a generalized chain of com- afunctionofhowmuchthepermutationdeviatesfrom ponentswheretheinteractionsbetweensitesorthecost theidentitypermutation. foranincorrectlyfilledsiteisnotconstantbutisdrawn In this model, we found that due to a state space fromadistributionofvalues. Suchasystemofquenched which could not be factorized in a basis defined by disorderwouldcharacterizeapermutationglasswhich lattice sites, even the superficially non-interacting sys- maycontaininterestingresultsduetotheuniquenature tem can exhibit phase-like transitions, i.e., temperature ofthestatespace. dependent changes in the value of the order parame- Supposingitispossibletodefinemoreinterestingin- ter which do note exhibit the properties typically as- teractions models, a natural investigation would con- sociated which phase transitions in infinite systems. cerntherenormalizationgrouppropertiesofthesystem. When interactions are introduced through a quadratic Specifically, we would be interested in how would one meanfieldterm,thesystemiscapableofexhibitingfive sumoverspecificstates(ascharacteristicofarenormal- regimesofthermalbehavior,andischaracterizedbytwo izationgrouptransformation)whenthestatespaceofa transition-temperaturescorrespondingtovariousquasi- systemlookslike, phasetransitions. The introduced model provides us with a basic ex- N actlysolublesystemforcertaininteractionassumptions = (42) system i S ⊗S andthusprovidesaconcretemodel-basedunderstand- i=1 Y ingof asystemwith anon-factorizablestatespace. Be- i.e.,isnotfactorizablealonglatticesites. cause of its utility and the type of results obtained, the modeldeservestobesubjecttothestandardextensions Finally,toconnectthismodelofpermutationstoprob- 9 lemsmorerelevanttoproteindesignitwouldprovenec- AppendixA:AlternativeDerivationsofEq.(31) essarytoincorporatethepossibilityofrepeatedcompo- nentsorthebackgroundgeometryofalatticechain. 1. Hubbard-StratonovichApproach Were-deriveEq.(31)usingtheHubbard-Stratonovich method. Westartthis derivationassumingλ = λ ; 2 2 −| | wewilllaterseeourresultingfreeenergycanbeanalyt- icallycontinuedtotheλ >0case. 2 Forλ = λ ,thepartitionfunctionis 2 2 −| | N Z (β;λ ,λ )= g (j)e−βλ1j+β|λ2|j2/2N. (A1) N 1 2 N ACKNOWLEDGMENTS Xj=0 Then,applyingtheidentity MW thanks Verena Kaynig-Fittkau for her unpub- eβ|λ2|j2/2N = N ∞ dxe−Nx2/2β|λ2|−jx, (A2) lGisihlseodn,wAorbkigtahilatPilnusmpmireedr, tVhiisnoitnhvaenstiMgaatnioonhaarnadn,Aamndy s2πβ|λ2|Z−∞ MichaelBrennerforhelpfuldiscussions. wehave Z (β;λ ,λ )= N ∞ dxe−Nx2/2β|λ2| N g (j)e−j(βλ1+x) N 1 2 N s2πβ|λ2|Z−∞ j=0 X = N ∞ dxe−Nx2/2β|λ2|Z (βλ +x), (A3) N 1 s2πβ|λ2|Z−∞ whereZ (x) N g (j)e−jx. FromEq.(11)wefound N ≡ j=0 N P ∞ Z (x)= dse−s 1+(s 1)e−x N. (A4) N − Z0 (cid:0) (cid:1) SoEq.(A3)becomes Z (β;λ ,λ )= N ∞ds ∞ dx 1+(s 1)e−βλ1−x Ne−s−Nx2/2β|λ2|. (A5) N 1 2 s2πβ|λ2|Z0 Z−∞ − (cid:0) (cid:1) Thefunctiontowhichweapplysteepestdescentisthen Solvingforsinthefirstequationwehave Nx2 s=N +1 eβλ1+x (A9) h(s,x)=s+ Nln 1+(s 1)e−βλ1−x . (A6) − 2β λ − − 2 | | andwiththesecondequationwehavethecondition (cid:0) (cid:1) Computingtheconditionsfor∂ h(s=s,x=x)=0and s x 1 ∂xh(s=s,x=x)=0weobtain,respectively, = (s 1). (A10) β λ −N − 2 | | e−βλ1−x 1 N =0 (A7) Substitutingthesecondconditionintothefirstyields − 1+(s 1)e−βλ1−x − x (s 1)e−βλ1−x s 1=N eβλ1−β|λ2|(s−1)/N, (A11) + − =0. (A8) β λ2 1+(s 1)e−βλ1−x − − | | − 10 or Wecanalsodefine e−β|λ2|(s−1)/N =eβλ1(N (s 1)). (A12) N − − (j) = (j)e−βλ0j, (A20) 0,N hO i O Thus,thesolutionforscanbeexpressedintermsofthe j=0 X LambertWfunctionas as the average with respect to our variational Hamilto- s 1 1 β λ nianEq.(A18). Thus,Eq.(A16)becomes − =1+ W | 2|eβλ1−β|λ2| . (A13) N β λ − N | 2| (cid:18) (cid:19) 1 λ F[ ] lnZ (βλ )+(λ λ ) j + 2 j2 N 0 1 0 0,N 0,N For λ > 0, we can employ the complex version of the H ≤−β − h i 2Nh i 2 Hubbard-Stratonovichidentity: f(λ ) (A21) 0 ≡ e−βλ2j2/2N = N ∞ dxe−Nx2/2βλ2−ijx. (A14) pDuiffteertehnetiamtianxgimfuwmitohfrtehsipseqctuatontλit0y.allGoiwvesnusjto com=- s2πβλ2 Z−∞ ∂lnZ (βλ )/∂(βλ ),wethenfind h i0,N N 0 0 − Working through an analogous steepestdescent proce- ∂ λ ∂ dure,wefindthattheequilibriumvalueforsis f′(λ )=(λ λ ) j + 2 j2 , (A22) 0 1 0 0,N 0,N − ∂λ h i 2N ∂λ h i 0 0 s 1 1 βλ j − =1 W 2eβλ1+βλ2 , (A15) whichifwetaketobezeroatsomeλ0 =λ0givesus N − βλ N ≡ N 2 (cid:18) (cid:19) ∂ λ ∂ which could have been extrapolated from Eq.(A13) by 0=(λ1 λ0) j 0,N + 2 j2 0,N . ltiabkriinugm|λc2o|n→diti−oλn,2.aFnrdomasthains eaxnparloesgsyiowniftohr tthheeenqouni-- − ∂λ0h i (cid:12)(cid:12)λ0=λ0 2N ∂λ0h i (cid:12)(cid:12)λ0(=Aλ203) To compute these deriv(cid:12)atives we make use of (cid:12)various interactingcase,itturnsouttheorderparameterinthis identities. Firstwenote caseisnotsbutrathers 1. − 1 ∂2 j2 = Z (βλ ), (A24) h i0,N Z (βλ )∂(βλ )2 N 0 N 0 0 so ∂ ∂2 j = lnZ (βλ ) 2. Gibbs-BogoliubovInequalityDerivationofEq.(31) ∂(βλ )h i0,N −∂(βλ )2 N 0 0 0 1 ∂2 = Z (βλ ) −Z (βλ )∂(βλ )2 N 0 We re-derive Eq.(31) using the Gibbs-Bogoliubov N 0 0 Inequality[10]. Theinequalityis 1 ∂ 2 + Z (βλ ) Z (βλ )2 ∂(βλ ) N 0 F[ ] F[ ]+ . (A16) N 0 (cid:18) 0 (cid:19) H ≤ H0 hH−H0i0 = j2 + j 2 . (A25) −h i0,N h i0,N TheHamiltonianwhichdefinesoursystemis Thislastequalityimplies N λ λ H=λ1Xi=1Iθi6=ωi + 2N2 Xi,j Iθi6=ωiIθj6=ωj ≡λ1j+ 2N2 j2, ∂(β∂λ0)hj2i0,N =−∂∂2(hβjλi00,)N2 +2hji0,N∂∂h(jβiλ0,0N), (A26) (A17) andourvariationalHamiltonianisinstead andsoEq.(A23)becomes N λ ∂ 2 H0 =λ0i=1Iθi6=ωi =λ0j. (A18) 0=(cid:20)λ1−λ0+ Nhji0,N(cid:21)∂λ0hji0,N(cid:12)λ0=λ0 X λ ∂2 j (cid:12) 2 h i0,N (cid:12) (A27) FromEq.(11)weknow − 2Nβ ∂λ20 λ0=λ0 (cid:12) (cid:12) 1 Tocomputethesequantities(cid:12)weneedtoapproximatethe F[ ]= lnZ (βλ ) 0 N 0 H −β partitionfunctionforourvariationalsystem. Usingthe = 1 ln ∞dse−s 1+(s 1)e−βλ0 N . −β − (cid:26)Z0 (cid:27) (cid:0) (cid:1) (A19)