Evolutionary Dynamics in Gene Networks and Inference Algorithms

Compu a ion 2015, 3, 99-113; doi:10.3390/compu a ion3010099 compu a ion ISSN 2079-3197 www.mdpi.com/jou nal/compu a ion Re iew E olu iona y Dynamics in Gene Ne wo ks and In e ence Algo i hms Daniel Aguila -Hidalgo 1,*, Ma ía C. Lemos 2 and An onio Có doba 2 1 Max Planck Ins i u e o he Physics o Complex Sys ems, Nö hni ze S aße 38, 01187 D esden, Ge many 2 Depa amen o de Física de la Ma e ia Condensada, Uni e sidad de Se illa, 41012 Se illa, Spain; E-Mails: [email p o ec ed] (M.C.L.); [email p o ec ed] (A.C.) * Au ho o whom co espondence should be add essed; E-Mail: daguila @pks.mpg.de; Tel.: +49-351-871-2120. Academic Edi o s: Ma nix Medema and Raine B ei ling Recei ed: 31 Oc obe 2014 / Accep ed: 3 Ma ch 2015 / Published: 13 Ma ch 2015 Abs ac : Dynamical in e ac ions among se s o genes (and hei p oduc s) egula e de elopmen al p ocesses and some dynamical diseases, like cance . Gene egula o y ne wo ks (GRNs) a e di ec ed ne wo ks ha de ine in e ac ions (links) among di e en genes/p o eins in ol ed in such p ocesses. Gene ic egula ion can be modi ied du ing he ime cou se o he p ocess, which may imply changes in he nodes ac i i y ha leads he sys em om a speci ic s a e o a di e en one a a la e ime (dynamics). How he GRN modi ies i s opology, o p ope ly d i e a de elopmen al p ocess, and how his egula ion was acqui ed ac oss e olu ion a e ques ions ha he e olu iona y dynamics o gene ne wo ks ackles. In he p esen wo k we e iew impo an me hodology in he ield and highligh he combina ion o hese me hods wi h e olu iona y algo i hms. In ecen yea s, his combina ion has become a powe ul ool o i models wi h he inc easingly a ailable expe imen al da a. Keywo ds: e olu iona y dynamics; e olu iona y algo i hms; gene egula o y ne wo ks OPEN ACCESS Compu a ion 2015, 3 100 1. In oduc ion Du ing e olu ion, o ganisms adap o he en i onmen o su i al ollowing na u al selec ion p ocesses, which in ol e modi ica ion o physiological cha ac e is ics o e gene a ions. Those modi ica ions imply changes in he genome ha may lead o di e en gene ic in e ac ions o pe o m a speci ic unc ion, i.e., gene egula ion. An e olu iona y p ocess can be di ided in wo pa s. The i s pa is he modi ica ion o gene egula ion ( ecombina ion ules and mu a ion), which leads o he modi ica ion o ce ain cha ac e is ics in a speci ic way. The second pa op imizes he cha ac e is ics modi ica ion by su i al o he bes -adap ed cha ac e is ics o hei pu pose (selec ion and i ness). F om a dynamical sys ems pe spec i e, hese p ocesses de ine an e olu iona y dynamics ha d i es a sys em om an ini ial s a e, a ac o A, which is s able unde ce ain en i onmen al condi ions, h ough di e en ajec o ies in a s a e space un il i eaches a new a ac o B, which is s able unde new en i onmen al condi ions (Figu e 1) [1]. Figu e 1. Scheme o e olu iona y dynamics om h ee di e en pe spec i es. (Top) F om he pe spec i e o dynamical sys ems; (Middle) F om he pe spec i e o he e olu iona y p ocess; (Bo om) F om he pe spec i e o ne wo ks opology. Compu a ion 2015, 3 101 In Figu e 1 op, he sys em in a ac o A modi ies i s s a e acco ding o ce ain dynamics ha may lead he sys em o di e en a ac o s. These a ge a ac o s may ep esen s a es, which a e a om he op imum possible s a e (in g ey). A selec ion owa ds an op imum may lead a ac o A o he op imum (o close o he op imum) a ac o B. Dashed ci cles ep esen he obus ness icini y whe e small pe u ba ions main ain he sys em unde he same a ac o (canaliza ion). Dynamical ules and s ochas ic pe u ba ions de ine he ajec o ies in he s a e space and de e mine hei inal s a es. The s a e space in he dynamical sys em ma ches wi h he i ness space in he e olu iona y p ocess. In Figu e 1middle, he sys em in A (g een do ) e ol es acco ding o ce ain dynamics and i can each a new s able s a e wi h low i ness (local maxima ep esen ed wi h g ey do s). A selec ion p ocess owa ds an op imum may lead A o a global maximum- i ness, o su icien ly high- i ness, B (blue do ). E olu iona y p ocesses inhe en ly de ine an op imiza ion mechanism so ha he app oxima ion o a ac o B is no a andom ajec o y. Heu is ic op imiza ion me hods and e olu iona y algo i hms mimic na u al selec ion using ecombina ion and mu a ion ules o app oxima e o a ac o B. The applica ion o heu is ic op imiza ion o he o mal de ini ion o he e olu iona y dynamics o gene ic egula ion may help o unde s and he dynamical ules ha lead om A o B. In p ac ice, e olu iona y dynamics o malisms use a se o gi en ules o e eal he inal s a es o a modi ied gene ic egula ion. In he ecen yea s, his app oach has been especially s udied in he ixa ion p obabili y o dele e ious mu a ions in i al quasispecies [2,3], and cance de elopmen [4–6]. The applica ion o e olu iona y algo i hms o he e olu iona y dynamics o gene ic egula ion de e mines he modi ica ion ules gi en bo h ini ial A and inal B s a es. In a p e ious wo k [7], we ha e ollowed his app oach o de e mine possible dynamical ules, which change gene ic egula ion be ween wo di e en de elopmen al s ages o a single de elopmen al p ocess [7]; Ge s ung e al. iden i ied possible dynamical ules o cance p og ession be ween di e en disease s ages [8]. In he p esen wo k, we e iew he ecen me hodology o e olu iona y dynamics and heu is ic op imiza ion, and discuss he applicabili y o such me hods o eal biological p oblems. Basic Concep s Gene egula ion o malisms. Gene egula ion includes a wide ange o mechanisms ha de e mine he concen a ion o gene p oduc s in a speci ic p ocess ha can be associa ed wi h he esul o a speci ic pheno ype o biological unc ion. In a simplis ic de ini ion, gene egula ion can be ep esen ed by he in e ac ion among ce ain genes (o gene p oduc s), which con ibu e o de elop biological p ocesses. Unde his de ini ion gene ic egula ion can be s udied in ne wo ked amewo ks, whe e he nodes ep esen genes (o hei p oduc s) and he links ep esen he egula o y in e ac ions among genes. In gene egula o y ne wo ks (GRN) we can ma ch nodes wi h exp essed (non-exp essed) genes ha a e unc ionally ac i e (inac i e) in a speci ic p ocess. Links can de ine posi i e egula o y in e ac ions, hey a o he ac i a ion o a non-exp essed gene o he ac i i y o an al eady exp essed one, o nega i e, hey a o he deac i a ion o an exp essed gene o block he ac i a ion o a non-exp essed one. The simples and mos common ep esen a ion o GRNs is he Boolean ype [7,9–33], whe e wo disc e e alues de ine he s a e o ne wo k nodes; common alues a e 1 and 0 (o +1 and −1) whe he he s a e is ac i e and inac i e espec i ely. Ex ensions o his ep esen a ion Compu a ion 2015, 3 102 comp ise mo e han wo disc e e alues o de ine g adual ac i i y o a node; his is usually conside ed in combina ion wi h Fuzzy-logic echniques [34]. Links ac ion can also be de ined as Boolean so ha alue 1 indica es a link be ween wo nodes, and 0 indica es ha hose wo nodes a e unconnec ed. The GRN link- alues de ine he adjacency ma ix M, whe e he en y mij desc ibes he in e ac ion be ween genes i and j; he diagonal en y mii desc ibes he sel -in e ac ion o gene i. I is also common o designa e an ac ion-in ensi y (weigh ) o he link using a g aded- alue in e al, hus he en y wij in he weigh ma ix W de ines he s eng h o he in e ac ion be ween gene i and j. A GRN is a di ec ed ne wo k, so M and W a e non-symme ic. Usually, GRNs a e spa se hough hey may con ain highly connec ed nodes (hubs). Con inuum o malisms p o ide mo e de ail ep esen a ions o GRNs. As a ecen ly explo ed ins ance, S-Sys ems de ines he nodes in a GRN as a iables in an o dina y di e en ial equa ion (ODE) sys em [13,35–42], which dynamics is de ined by: 11 ij ij nn g h i ijij jj dX X X d == =α −β ∏∏ (1) whe e Xi de ine he genes exp ession (ne wo k nodes), which is gi en in ime as a unc ion o he exp ession le els o he o he genes. The dynamical eac ion kine ics in Equa ion (1) is de ined using a powe -law o malism ha allows cap u ing complex non-linea ela ions. The wo eac ion e ms in Equa ion (1) o mally desc ibe syn hesis and deg ada ion o gene i, which a e in luenced by genes j. αi and βi a e eac ion a es (syn hesis and deg ada ion espec i ely) and gij and hij a e he posi i e kine ic o de s o he eac ions ha indica e he s eng h o he in luence o gene j in he syn hesis and deg ada ion o gene i. O he ins ances o GRN de ini ion a e based on linea ODE sys ems [16,43–48], a i icial neu al ne wo ks (ANN) [29,34,49] and e icula sys ems [30–32]. E olu iona y dynamics ules. The opology o a GRN de ines he exp ession o he genes (nodes ac i i y) in a de e mined ime poin . E olu iona y ules can modi y nodes ac i i y and ne wo k wi ing in ime so ha he inal opology o he ne wo k may be di e en om he ini ial one. In disc e e sys ems, e olu ion ules a e gene ally de ined as h eshold- ype unc ions, ei he s ep-like (Hea iside) o sigmoidal (like he Hill unc ion); he ules a e applied using di e ence equa ions wi h a disc e e ime coun e . In con inuum ep esen a ions he ules a e applied using di e en ial equa ions and a eal a iable o he ime. Du ing e olu ion, GRNs can be anked acco ding o a p ede ined c i e ion using a i ness unc ion. This unc ion measu es he adap abili y o he e ol ed sys em and d i es he e olu ion p ocess o an a ac o . The i ness unc ion is de ined acco ding o he speci ic sys em o s udy. In Figu e 1bo om, Ne wo k A (wi h g een nodes) ep esen s a ce ain unc ion o a sys em (pheno ypic/geno ypic cha ac e is ic in e ms o gene ne wo ks). Following a speci ic dynamics, ne wo k A changes i s opology by ewi ing he links and/o changing he numbe o nodes. This new opology may ep esen a s a e a om he op imum (ne wo ks wi h nodes in g ey). A selec ion p ocess owa ds an op imum may lead ne wo k A o ne wo k B (wi h blue nodes) sa is ying a speci ic unc ionali y. Compu a ion 2015, 3 103 2. E olu iona y Dynamics Me hods Ou s a ing poin in he e olu iona y dynamics me hods is Wagne ’s model [9,10]. Wagne s udied he in luence o gene duplici y on gene exp ession in me azoan de elopmen . The GRN is disc e ely ep esen ed using alue 1 o ac i e nodes and −1 o inac i e nodes. In Equa ion (2), Si( ) ep esen s he s a e o he nodes in ime . The sys em e ol es applying he di e ence equa ions: () () 1 N iijj j S wS =  +τ =    (2) whe e (x) indica es he sign o he exp ession, and 1 () N ij j j x wS = =. The exp ession o he gene Si in ime ( + τ) is a unc ion o he ac ion o he genes Sj in ime wi h wij he egula o y s eng h o he in e ac ion o he gene Sj on Si. I x < 0 hen (x) = −1 and he egula ion is nega i e, hus he a ge gene is deac i a ed. I x > 0 hen (x) = +1 and he egula ion is posi i e, hus he a ge gene is ac i a ed. I x = 0 hen (x) = 0 and he a ge node is no modi ied. The cha ac e is ic ime τ indica es he ime uni du ing he e olu iona y p ocess. This dynamics leads ei he o an equilib ium o oscilla o y s a e. Fo simplici y, Wagne s udied only he equilib ium s a e ˆ S and analyzed he epigene ic s abili y ha is he numbe o s able s a es. An op imal s a e Sop exis s du ing he de elopmen al p ocess. I he inal equilib ium s a e is di e en om Sop , he de elopmen al p ocess may su e modi ica ions and, he e o e, he i ness o he adul o ganism is educed. The e olu iona y p ocess, s a ing om an ini ial s a e, is nume ically simula ed o a se o indi iduals and is cha ac e ized by a ce ain s a e S. The dis ance o he indi idual i ness alue agains he op imal is measu ed using he Hamming dis ance: 1 11 ˆ [,] 22 N op ii i dSS SS N= =−  (3) F om Exp ession (3), he i ness is gi en by he unc ion: 2 ˆ [,] exp 2 dSS F =−    (4) whe e β is a posi i e pa ame e ha measu es he selec ion s eng h. To es whe he he e olu ion ule leads o an equilib ium s a e, i is applied wice. I he ne wo k is in an equilib ium s a e he i ness is measu ed, o he wise i is assigned he minimum i ness alue, exp(−1/β). The indi iduals measu ed a e anked by hei i ness alue. The ecombina ion p ocess is only pe o med o e selec ed indi iduals. Two indi iduals a e selec ed o ecombina ion i hei no malized i ness is highe han a andom numbe wi h uni o m dis ibu ion be ween 0 and 1. The non-selec ed indi iduals a e disca ded. The ecombina ion p ocess implies he exchange o ows in he weigh ma ix W be ween wo indi iduals selec ed acco ding o a highes i ness ule. The weigh s wij ≠ 0 ha e a p obabili y o modi ica ion (mu a ion) based on a pseudo andom alue. Wagne de e mined ha gene duplici y does no a ec he exp ession le el. Compu a ion 2015, 3 104 Many au ho s ha e conside ed a ia ions on Wagne ’s model. As ins ances, Siegal and Be gman p oposed he same dynamical law as in Equa ion (2). As a di e ence, he au ho s in oduce a sign unc ion in he o m: 2 () 1 1ax x e− =− + (5) The sign unc ion (x) is a sigmoidal unc ion wi h a con olling he speed o change om an ac i a ion o a ep ession s a e [15]. Wi h his model, he au ho s e isi ed Wadding on’s canaliza ion concep ha is he obus ness du ing de elopmen o changes in he genome [50]. They claim ha he dynamical GRN d i ing a de elopmen al p ocess, desc ibed by i e a ion o Equa ions (2) and (5), e ol es o cons ain he gene ic sys em o p oduce canaliza ion, e en wi hou a selec ion owa ds an op imum i ness [15]. To de ine he inal s a e, ˆ S, o he e olu iona y p ocess, hey de ine a measu e simila o a a iance in he o m: () () τ 1 () ( ), () S DS S θ= − ψ= θ τ (6) whe e () () 2 1 1 (), () 4 N ii i DS S s s N= θ= −  (7) and S is he mean exp ession le el du ing he ime in e al ( − τ, …, ). The ini ial s a e is andomly selec ed and he s eady s a e, ˆ S, is eached when ψ(ˆ S) < 10−4. The sys em is conside ed uns able i his h eshold is no each wi hin 100 i e a ions. In his model, he i ness o an uns able indi idual is 0. The i ness unc ion in s eady-s a e sys ems is de ined as: () () ˆ, ˆexp op DSS FS     =−   β   (8) wi h β he selec ion s eng h and Sop p ede ined he op imum solu ion. No e ha he dynamical law, Equa ion (2), and Hamming dis ance, Equa ion (3), de ined in Wagne ’s model is a special case o his one conside ing he limi o a∞. The ecombina ion is pe o med by andomly exchanging W ows in indi iduals wi h simila i ness. Mu a ions a e pe o med acco ding o a no mal dis ibu ion. The new indi idual is conside ed pa o he new popula ion i i eaches he s eady s a e and i s i ness is highe han a andom numbe wi h uni o m dis ibu ion. O he ins ances s udy he e ec o mu a ions wi h espec o he s abili y o he e ol ed s a e. Masel s udied gene ic assimila ions ha a e mu a ions p e iously induced by ex e nal en i onmen al condi ions ha la e u n in o inhe i ed gene ic modi ica ions [17]. The au ho added a noise componen ε in he o m o a andom numbe d own om a no mal dis ibu ion. The noise componen may al e he s a e o a node si in he dynamic equa ion: () ( ) () 1 ii s WS =−+ε (9) Compu a ion 2015, 3 105 whe e (x) = 0 i x < 0 and (x) = 1 o he wise wi h x = WS( − 1)i + ε. I S does no emain a some s able alue du ing ou ime s eps in a ow wi hin 100 ime s eps, hen he indi idual, ep esen ed by he ma ix W, does no each equilib ium and i is assumed o be un iable (W is disca ded). The au ho ound ha gene ic assimila ion could occu in absence o a selec ion p ocess ha a o s he assimila ing modi ica ion. This is also in ag eemen wi h Wadding on’s canaliza ion. Simila ly, Espinosa-So o e al. analyzed he e ec o gene ic and non-gene ic pe u ba ions on an indi idual and i s ela ion wi h he pheno ypic plas ici y ha a egula o y ci cui may each o ease he adap i e e olu ion [31]. Aze edo e al. de eloped a model in which sexual ep oduc ion p oduces nega i e epis asis ha inc eases he capaci y o pu ge deadly mu a ions. This sugges s ha sexual ep oduc ion selec s condi ions o a o i s own main enance [21]. Kaneko e al. s udied he obus ness agains noise and mu a ions in gene ic exp ession. They used a di e en ype o model based on s ochas ic di e en ial equa ions [44]. 3. In e ence o GRNs Using E olu iona y Algo i hms One o he mo i a ions o e olu iona y dynamics me hods is he in e p e a ion o eal expe imen al da a. Se e al me hods a e used o in e his in e p e a ion. In he case o genomics, his in e ence p o ides GRNs ha a e able o sa is y he expe imen al condi ions. The mos ecen ly used me hods a e e olu iona y algo i hms (EA) [51]. Among he EA one can ind he gene ic algo i hms (GA), e olu ion s a egies (ES), e olu iona y p og amming (EP), simula ed annealing (SA), an colonies (AC) and immune algo i hms (IA) [51]. Mo eo e , hyb id algo i hms a e also conside ed. He e, wo o mo e me hods a e combined, o hey a e applied wi h o he selec ion/lea ning me hods such as a i icial neu al ne wo ks (ANN) and S-sys ems. We i s discuss wo o he mos common single me hods. One o he, so-called, classical e olu iona y algo i hms a e he gene ic algo i hms. The GAs a e based on na u al selec ion. They con empla e inhe i ance ( ecombina ion), mu a ion, adap a ion and eli ism. The gene al scheme (Figu e 2) o hese me hods is he ollowing: he s a ing poin is he gene a ion o indi iduals (o en called ch omosomes) ep esen ing he possible solu ions o a p oblem. The indi iduals a e anked acco ding o a speci ic i ness unc ion ha measu es he p oximi y o he indi idual’s solu ion o a p ede ined op imal solu ion (commonly he i o expe imen al da a). The selec ion o indi iduals o o m o sp ing ha a e new indi iduals, which a e supposed o main ain mos o he good quali ies o hei pa en s, is i ness dependen . The o sp ing is o m by c osso e me hods o he selec ed indi iduals. The new indi iduals can mu a e. Mu a ions a e small pe u ba ions, which may e en ually lead o be e cha ac e is ics. This ep oduc i e cycle in ends o imp o e he quali y o he indi iduals acco ding o he i ness c i e ia. In he g oup o he conside ed mode n EA, one o he mos used is he simula ed annealing. While GA mimic laws o e olu ion in biology, SA emula es a me allu gic echnique used o ob ain uni o m composi ion alloys. The i ness (ene gy) unc ion has mul iple local minima (me as able s a es) and i should be minimal in he he modynamic equilib ium. In expe imen al annealing, a solid is hea ed o a high empe a u e, whe eby pa icles in he liquid phase ha e a g ea eedom o adap o a minimum ene gy con igu a ion. Then, he cooling p og am slowly dec eases empe a u e T, so ha o each T alue he sys em eaches an equilib ium s a e. The p ocess can be desc ibed s a ing om he Compu a ion 2015, 3 106 maximum empe a u e and conside ing ha o each T alue he solid eaches i s he mal equilib ium. Each empe a u e-dependen s a e is cha ac e ized by an ene gy E gi en by he Bol zmann dis ibu ion, p o ided ha he cooling p og am is su icien ly slow. As empe a u e lowe s, he Bol zmann dis ibu ion concen a es in lowe ene gy s a es. Finally, when he empe a u e is close o ze o, only he lowes ene gy s a e has a nonze o p obabili y. I he cooling p og am is done oo as , he solid can each me as able s uc u es. Figu e 2. Wo k low schemes o gene al gene ic algo i hms (GA) and simula ed annealing (SA) implemen a ion. In heu is ic op imiza ion, he SA me hod ma ches he i ness unc ion wi h he ene gy unc ion, which includes an a i icial con ol pa ame e as empe a u e. The se o possible solu ions is mapped Compu a ion 2015, 3 107 o he possible s a es o he physical sys em. The s eady s a e ( he op imal solu ion) ha minimizes he ene gy ( i ness unc ion) is achie ed wi h high p obabili y a he end o he SA p ocess. The empe a u e pa ame e g adually dec eases om an ini ial high alue. Ins ances o hese me hods in e olu iona y dynamics o GRN a e he ollowing: Ina p e ious wo k we applied a GA o deciphe he dynamical ules ha d i e he e olu ion o a GRN be ween wo di e en de elopmen al s ages. The GA modi ies bo h he ne wo k links and he e olu ion ules applied o each node. We de e mined ha , in many cases, di e en ules migh p oduce he same e ec in he ne wo k e olu ion [7]. Kobayashi e al. applied a SA me hod o op imize a ne wo k opology so ha he GRN exhibi s able sus ained oscilla ions wi h a p ede ined pe iod [46,47]. The au ho s conside ed he Me opolis algo i hm om Mon e Ca lo me hod wi h he cooling p og am T = μF, whe e T is he e ec i e empe a u e, μ is a no maliza ion cons an and F is he i ness unc ion. The i ness F dec eases wi h T un il i eaches a su icien ly low alue. The i ness unc ion is de ined as: () 22 0 22 0 PP FPP −σ =+ (10) whe e P0 is he sea ched pe iod, P is he mean pe iod du ing he simula ion ime and σ is he pe iod P a iance, which should anish i he GRN dynamics is eally pe iodic. Though single me hods can ge su icien ly good solu ions, mos ecen s udies a e based on hyb id me hods. Some o he mos ele an a e discussed below. Tominaga e al. designed an in e ence me hod o GRNs combining S-sys ems wi h GA [35]. The mo i a ion o his esea ch is o i a dynamical GRN ou pu wi h expe imen al ime se ies. In Tominaga’s me hod, he se o S-sys em pa ame e s (αi, βi, gij, hij) o ms an indi idual. The i ness unc ion (Equa ion (11)) measu es he p oximi y o he S-sys em solu ion o each indi idual wi h espec o an expe imen al empo al se ies. () () () 1 2 ' 11 nm ii ij i X X FX − ==   −  =      (11) whe e n is he numbe o obse able s a e a iables, m is he numbe o expe imen al poin s, ' i X ( ) a e he simula ed alues and Xi( ) a e he expe imen al alues. Kikuchi e al. pe o med modi ica ions on Tominaga’s model o inc ease bo h he numbe o pa ame e s o op imize and he con e gence a e o an op imum i , Equa ion (12). This modi ica ion adds an ex a e m ha a o s he indi iduals wi h low gij and hij. () () () 1 2 ' 11 , ,, nm ii ij ij i j ij iji j i X X Fcnmgh X − == ≠    −  =++           (12) whe e c is a balance weigh be ween he wo e ms [36]. The indi iduals wi h gij and hij unde a p e-de ined h eshold a e se o ze o. This in oduces a p uning mechanism ha makes he ne wo k successi ely spa se . A e a ce ain numbe o gene a ions minimizing he ne wo k opology, Equa ion (11) is used o a e inemen sea ch.