CG920 Genomics
Lesson 11
Systems Biology
Jan Hejátko
Functional Genomics and Proteomics of Plants,
CEITEC - Central European Institute of Technology
And
National Centre for Bimolecular Research,
Faculty of Science,
Masaryk University, Brno
hejatko@sci.muni.cz, www.ceitec.eu
2
2
 Literature sources for Chapter 12:
 Wilt, F.H., and Hake, S. (2004). Principles of Developmental Biology. (New York ; London: W. W.
Norton)
 Eden, E., Navon, R., Steinfeld, I., Lipson, D., and Yakhini, Z. (2009). GOrilla: a tool for discovery
and visualization of enriched GO terms in ranked gene lists. BMC Bioinformatics 10, 48.
 The Arabidopsis Genome Initiative. (2000). Analysis of the genome sequence of the flowering
plant Arabidopsis thaliana. Nature 408, 796-815.
 Benitez, M. and Hejatko, J. Dynamics of cell-fate determination and patterning in the vascular
bundles of Arabidopsis thaliana (submitted)
 de Luis Balaguer MA, Fisher AP, Clark NM, Fernandez-Espinosa MG, Moller BK, Weijers D,
Lohmann JU, Williams C, Lorenzo O, Sozzani R. 2017. Predicting gene regulatory networks by
combining spatial and temporal gene expression data in Arabidopsis root stem cells. Proc Natl
Acad Sci U S A 114(36): E7632-E7640.
Literature
3
3
 Definition of Systems Biology
 Tools
 Gene Ontology Analysis
 Bayesian Networks
 Molecular/Gene Regulatory Networks Modeling
 Inferring Gene Regulatory Networks from Large Omics
Datasets
Outline
4
4
Systems biology is the computational and mathematical
analysis and modeling of complex biological systems. It
is a biology-based interdisciplinary field of study that
focuses on complex interactions within biological
systems, using a holistic approach (holism instead of the
more traditional reductionism) to biological research
(Wikipedia).
Definition
5
5
Systems biology is the study of biological systems
whose behaviour cannot be reduced to the linear sum of
their parts’ functions. Systems biology does not
necessarily involve large numbers of components or vast
datasets, as in genomics or connectomics, but often
requires quantitative modelling methods borrowed from
physics (Nature).
Definition
6
6
Nice explanatory video by Dr. Nathan Price,
associate director of the Institute for Systems Biology at
https://www.youtube.com/watch?v=OrXRl_8UFHU.
Definition
7
7
 Definition of Systems Biology
 Tools
 Gene Ontology analysis
Outline
8
Results of –omics Studies vs
Biologically Relevant Conclusions
□ Results of –omics studies represent huge amount of data,
e.g. genes with differential expression. But how to get any
biologically relevant conclusions out of it?
gene locus sample_1 sample_2 status value_1 value_2 log2(fold_change) test_stat p_value q_value significant
AT1G07795 1:2414285-2414967 WT MT OK 0 1,1804 1.79769e+308
1.79769e+3
08 6.88885e-05
0,00039180
1 yes
HRS1 1:4556891-4558708 WT MT OK 0 0,696583 1.79769e+308
1.79769e+3
08 6.61994e-06
4.67708e-
05 yes
ATMLO14 1:9227472-9232296 WT MT OK 0 0,514609 1.79769e+308
1.79769e+3
08 9.74219e-05
0,00053505
5 yes
NRT1.6 1:9400663-9403789 WT MT OK 0 0,877865 1.79769e+308
1.79769e+3
08 3.2692e-08
3.50131e-
07 yes
AT1G27570 1:9575425-9582376 WT MT OK 0 2,0829 1.79769e+308
1.79769e+3
08 9.76039e-06 6.647e-05 yes
AT1G60095 1:22159735-22162419 WT MT OK 0 0,688588 1.79769e+308
1.79769e+3
08 9.95901e-08
9.84992e-
07 yes
AT1G03020 1:698206-698515 WT MT OK 0 1,78859 1.79769e+308
1.79769e+3
08 0,00913915 0,0277958 yes
AT1G13609 1:4662720-4663471 WT MT OK 0 3,55814 1.79769e+308
1.79769e+3
08 0,00021683 0,00108079 yes
AT1G21550 1:7553100-7553876 WT MT OK 0 0,562868 1.79769e+308
1.79769e+3
08 0,00115582 0,00471497 yes
AT1G22120 1:7806308-7809632 WT MT OK 0 0,617354 1.79769e+308
1.79769e+3
08 2.48392e-06
1.91089e-
05 yes
AT1G31370 1:11238297-11239363 WT MT OK 0 1,46254 1.79769e+308
1.79769e+3
08 4.83523e-05
0,00028514
3 yes
APUM10 1:13253397-13255570 WT MT OK 0 0,581031 1.79769e+308
1.79769e+3
08 7.87855e-06
5.46603e-
05 yes
AT1G48700 1:18010728-18012871 WT MT OK 0 0,556525 1.79769e+308
1.79769e+3
08 6.53917e-05
0,00037473
6 yes
AT1G59077 1:21746209-21833195 WT MT OK 0 138,886 1.79769e+308
1.79769e+3
08 0,00122789 0,00496816 yes
AT1G60050 1:22121549-22123702 WT MT OK 0 0,370087 1.79769e+308
1.79769e+3
08 0,00117953 0,0048001 yes
Ddii et al., unpublished
AT4G15242 4:8705786-8706997 WT MT OK 0,00930712 17,9056 10,9098 -4,405231.05673e-05 7.13983e-05 yes
AT5G33251 5:12499071-12500433 WT MT OK 0,0498375 52,2837 10,0349 -9,8119 0 0 yes
AT4G12520 4:7421055-7421738 WT MT OK 0,0195111 15,8516 9,66612 -3,900439.60217e-05 0,000528904 yes
AT1G60020 1:22100651-22105276 WT MT OK 0,0118377 7,18823 9,24611 -7,503826.19504e-14 1.4988e-12 yes
AT5G15360 5:4987235-4989182 WT MT OK 0,0988273 56,4834 9,1587 -10,4392 0 0 yes
Excample of an output of transcriptional profiling study using Illumina sequencing
performed in our lab. Shown is just a tiny fragment of the complete list,
copmprising about 7K genes revealing differential expression in the studied
mutant.
8
9
Plant Vascular Tissue
Development
□ Vascular tissue as a developmental model for GO analysis and MRN
modeling
Lehesranta etal., Trends in Plant Sci (2010)
10
WT hormonal mutant
Hormonal Control Over Vascular
Tissue Development
□ Plant Hormones Regulate Lignin Deposition in Plant Cell
Walls and Xylem Water Conductivity
WT mutant
lignified cell walls
Water Conductivity
WT hormonal mutants
10
11
WT hormonal mutant
Hormonal Control Over Vascular
Tissue Development
□ Transcriptional profiling via RNA sequencing
mRNA
Sequencing by Illumina and
number of transcripts determination
mRNA
cDNA cDNA
11
12
Results of –omics Studies vs
Biologically Relevant Conclusions
□ Transcriptional profiling yielded more then 9K differentially
regulated genes…
gene locus sample_1 sample_2 status value_1 value_2 log2(fold_change) test_stat p_value q_value significant
AT1G07795 1:2414285-2414967 WT MT OK 0 1,1804 1.79769e+308
1.79769e+3
08 6.88885e-05
0,00039180
1 yes
HRS1 1:4556891-4558708 WT MT OK 0 0,696583 1.79769e+308
1.79769e+3
08 6.61994e-06
4.67708e-
05 yes
ATMLO14 1:9227472-9232296 WT MT OK 0 0,514609 1.79769e+308
1.79769e+3
08 9.74219e-05
0,00053505
5 yes
NRT1.6 1:9400663-9403789 WT MT OK 0 0,877865 1.79769e+308
1.79769e+3
08 3.2692e-08
3.50131e-
07 yes
AT1G27570 1:9575425-9582376 WT MT OK 0 2,0829 1.79769e+308
1.79769e+3
08 9.76039e-06 6.647e-05 yes
AT1G60095 1:22159735-22162419 WT MT OK 0 0,688588 1.79769e+308
1.79769e+3
08 9.95901e-08
9.84992e-
07 yes
AT1G03020 1:698206-698515 WT MT OK 0 1,78859 1.79769e+308
1.79769e+3
08 0,00913915 0,0277958 yes
AT1G13609 1:4662720-4663471 WT MT OK 0 3,55814 1.79769e+308
1.79769e+3
08 0,00021683 0,00108079 yes
AT1G21550 1:7553100-7553876 WT MT OK 0 0,562868 1.79769e+308
1.79769e+3
08 0,00115582 0,00471497 yes
AT1G22120 1:7806308-7809632 WT MT OK 0 0,617354 1.79769e+308
1.79769e+3
08 2.48392e-06
1.91089e-
05 yes
AT1G31370 1:11238297-11239363 WT MT OK 0 1,46254 1.79769e+308
1.79769e+3
08 4.83523e-05
0,00028514
3 yes
APUM10 1:13253397-13255570 WT MT OK 0 0,581031 1.79769e+308
1.79769e+3
08 7.87855e-06
5.46603e-
05 yes
AT1G48700 1:18010728-18012871 WT MT OK 0 0,556525 1.79769e+308
1.79769e+3
08 6.53917e-05
0,00037473
6 yes
AT1G59077 1:21746209-21833195 WT MT OK 0 138,886 1.79769e+308
1.79769e+3
08 0,00122789 0,00496816 yes
AT1G60050 1:22121549-22123702 WT MT OK 0 0,370087 1.79769e+308
1.79769e+3
08 0,00117953 0,0048001 yes
Ddii et al., unpublished
AT4G15242 4:8705786-8706997 WT MT OK 0,00930712 17,9056 10,9098 -4,405231.05673e-05 7.13983e-05 yes
AT5G33251 5:12499071-12500433 WT MT OK 0,0498375 52,2837 10,0349 -9,8119 0 0 yes
AT4G12520 4:7421055-7421738 WT MT OK 0,0195111 15,8516 9,66612 -3,900439.60217e-05 0,000528904 yes
AT1G60020 1:22100651-22105276 WT MT OK 0,0118377 7,18823 9,24611 -7,503826.19504e-14 1.4988e-12 yes
AT5G15360 5:4987235-4989182 WT MT OK 0,0988273 56,4834 9,1587 -10,4392 0 0 yes
Excample of an output of transcriptional profiling study using Illumina sequencing
performed in our lab. Shown is just a tiny fragment of the complete list,
copmprising about 7K genes revealing differential expression in the studied
mutant.
12
13
Gene Ontology Analysis
□ One of the possible approaches is to study gene ontology, i.e.
previously demonstrated association of genes to biological
processes
Ddii et al., unpublished
14
Gene Ontology Analysis
□ Several tools allow statistical evaluation of enrichment for
genes associated with specific processes
Eden et al., BMC Biinformatics (2009)
One of such recent and very useful tools is Gorilla software, freely available at
http://cbl-gorilla.cs.technion.ac.il/.
14
15
Gene Ontology Analysis
□ Several tools allow statistical evaluation of enrichment for
genes associated with specific processes
16
Gene Ontology Analysis
□ Several tools allow statistical evaluation of enrichment for
genes associated with specific processes
17
Gene Ontology Analysis
□ Several tools allow statistical evaluation of enrichment for
genes associated with specific processes
18
Gene Ontology Analysis
□ Several tools allow statistical evaluation of enrichment for
genes associated with specific processes
19
Gene Ontology Analysis
□ Several tools allow statistical evaluation of enrichment for
genes associated with specific processes
20
20
 Definition of Systems Biology
 Tools
 Gene Ontology analysis
 Bayesian Networks
Outline
21
Bayesian Networks
 What are Bayesian networks?
 Probabilistic Graphical Model that can be used to build models from
data and/or expert opinion
https://www.youtube.com/watch?v=4fcqyzVJwH
M
21
22
Bayesian Networks
 What are Bayesian Networks?
 Probabilistic Graphical Model that can be used to build models from
data and/or expert opinion
 can be used for a wide range of tasks including prediction, anomaly
detection, diagnostics, automated insight, reasoning, time series
prediction and decision making under uncertainty
 NODES
 each node represents a variable such as someone's height, age or gender.
A variable might be discrete, such as Gender = {Female, Male} or might be
continuous such as someone's age
 LINKS
 added between nodes to indicate that one node directly influences the
other
22
23
Bayesian Networks
NODES
LINKS/EDGES
23
24
Asia Bayesian Network
https://www.bayesserver.com/
24
25
25
 Definition of Systems Biology
 Tools
 Gene Ontology analysis
 Bayesian Networks
 Molecular/Gene Regulatory Networks Modeling
Outline
26
□ Vascular tissue as a developmental model for MRN modeling
Benitez and Hejatko, PLoS One, 2013
Molecular Regulatory
Networks Modeling
27
Molecular Regulatory Networks Modeling
□ Literature search for published data and creating small
database
Interaction Evidence References
A-ARRs –| CK
signaling
Double and higher order type-A ARR mutants show
increased sensitivity to CK.
Spatial patterns of A-type ARR gene expression and
CK response are consistent with partially redundant
function of these genes in CK signaling.
A-type ARRs decreases B-type ARR6-LUC.
Note: In certain contexts, however, some A-ARRs
appear to have effects antagonistic to other A-
ARRs.
[27]
[27]
[13]
[27]
AHP6 –| AHP ahp6 partially recovers the mutant phenotype of the
CK receptor WOL.
Using an in vitro phosphotransfer system, it was
shown that, unlike the AHPs, native AHP6 was
unable to accept a phosphoryl group. Nevertheless,
AHP6 is able to inhibit phosphotransfer from other
AHPs to ARRs.
[9]
[9]
Benitez and Hejatko, PLoS One, 2013
28
Molecular Regulatory Networks Modeling
□ Formulating logical rules defining the model dynamics
Network node Dynamical rule
CK 2 If ipt=1 and ckx=0
1 If ipt=1 and ckx=1
0 else
CKX 1 If barr>0 or arf=2
0 else
AHKs ahk=ck
AHPs 2 If ahk=2 and ahp6=0 and aarr=0
1 If ahk=2 and (ahp6+aarr<2)
1 If ahk=1 and ahp6<1
0 else
B-Type ARRs 1 If ahp>0
0 else
A-Type ARRs 1 If arf<2 and ahp>0
0 else
Benitez and Hejatko, PLoS One, 2013
29
Molecular Regulatory Networks Modeling
□ Specifying mobile elements and their model behaviour
According to experimental evidence for the system under study, the hormone IAA, the
peptide TDIF, and the microRNA MIR165/6 are able to move among the cells. In the
case of TDIF and MIR165/6, the mobility is defined as diffusion and is given by the
following equation:
g(t+1)T[i]= H(g(t)[i]+ D (g(t)[i+1]+g(t)[i-1] – N(g(t)[i]))-b) (2),
where g(t)T[i] is the total amount of TDIF or MIR165 in cell (i). D is a parameter that
determines the proportion of g that can move from any cell to neighboring ones and is
correlated to the diffusion rate of g. b is a constant corresponding to a degradation term.
H is a step function that converts the continuous values of g into a discrete variable that
may attain values of 0, 1 or 2. N stands for the number of neighbors in each cell.
Boundary conditions are zero-flux. In the case of IAA, the mobility is defined as active
transport dependent on the radial localization of the PIN efflux transporters and is
defined by the equation:
iaa(t+1)T[i]=Hiaa(iaa(t)[i]+Diaa(pin(t)[i+1])(iaa(t)[i+1])+Diaa(pin(t)[i-1])(iaa(t)[i-
1])−N(Diaa)(pin(t)[i])(iaa(t)[i])−biaa) (3),
where Diaa is a parameter that determines the proportion of IAA that can be transported
among cells. The transport depends on the presence of IAA and PIN in the cells and
biaa corresponds to a degradation term. As in equation 2, H is a step function that
converts the continuous values to discrete ones and N stands for the number of
neighbors in each cell. Boundary conditions for IAA motion are also zero-flux.
29
30
Molecular Regulatory Networks Modeling
□ Preparing the first version of the model and its testing
The proposed model considers data that we identified and evaluated through an
extensive search (up to January 2012). It takes into account molecular
interactions, hormonal and expression patterns, and cell-to-cell communication
processes that have been reported to affect vascular patterning in the bundles of
Arabidopsis. The model components and interactions are graphically presented in
the figure above. In the network model, nodes stand for molecular elements
regulating one another’s activities. Most of the nodes can take only 1 or 0 values
(light gray nodes in the figure), corresponding to “present” or “not present,”
respectively. Since the formation of gradients of hormones and diffusible
elements may have important consequences in pattern formation, mobile
elements TDIF and MIR, as well as members of the CK and IAA signaling
systems, can take 0, 1 or 2 values (dark gray nodes in the figure above) Benitez
and Hejatko, submitted.
30
31
Molecular Regulatory Networks Modeling
□ Specifying of missing interactions via informed predictions
Interaction Evidence References
A-ARRs –| CK
signaling
Double and higher order type-A ARR mutants show
increased sensitivity to CK.
Spatial patterns of A-type ARR gene expression and CK
response are consistent with partially redundant function of
these genes in CK signaling.
A-type ARRs decreases B-type ARR6-LUC.
Note: In certain contexts, however, some A-ARRs appear to
have effects antagonistic to other A-ARRs.
[27]
[27]
[13]
[27]
AHP6 –| AHP ahp6 partially recovers the mutant phenotype of the CK
receptor WOL.
Using an in vitro phosphotransfer system, it was shown that,
unlike the AHPs, native AHP6 was unable to accept a
phosphoryl group. Nevertheless, AHP6 is able to inhibit
phosphotransfer from other AHPs to ARRs.
[9]
[9]
CK → PIN7 radial
localization
Predicted interaction (could be direct or indirect)
Informed by the following data:
During the specification of root vascular cells in Arabidopsis
thaliana, CK regulates the radial localization of PIN7.
Expression of PIN7:GFP and PIN7::GUS is upregulated by
CK with no significant influence of ethylene.
In the root, CK signaling is required for the CK regulation of
PIN1, PIN3, and PIN7. Their expression is altered in wol,
cre1, ahk3 and ahp6 mutants.
[18]
[18,20]
[19]
CK→ APL Predicted interaction (could be direct or indirect)
Consistent with the fact that APL overexpression prevents or
delays xylem cell differentiation, as does CKs.
Partially supported by microarray data and phloem-specific
expression patterns of CK response factors.
[21]
(TAIR,
ExpressionSet:
1005823559,
[22])
32
Molecular Regulatory Networks Modeling
□ Preparing the next version of the model and its testing
Benitez and Hejatko, PLoS One, 2013
In comparison to the model shown on slide 21, the final version of the model
contains the predicted interactions (dashed lines).
32
33
□ Good model should be able to simulate reality
Benitez and Hejatko, PLoS One, 2013
Molecular Regulatory
Networks Modeling
34
□ Formulating equations describing the relationships in the model
Molecular Regulatory
Networks Modeling
Static nodes: gn(t+1)=Fn(gn1(t),gn2(t),..., gnk(t))
Mobile nodes: g(t+1)T [i]= H(g(t) [i]+ D (g(t) [i+1]+g(t) [i-1] – N(g(t) [i]))-b)
state in the time t+1
state in the time tlogical rule function
state in the time t+1
Amount if TDIF or MIR165 in cell i
proportion of movable element
constant corresponding to a degradation term
34
35
□ Good model should be able to simulate reality
Benitez and Hejatko, submitted
Molecular Regulatory
Networks Modeling
Static nodes: gn(t+1)=Fn(gn1(t),gn2(t),..., gnk(t))
Mobile nodes: g(t+1)T [i]= H(g(t) [i]+ D (g(t) [i+1]+g(t) [i-1] – N(g(t) [i]))-b)
The initial conditions specify the initial state of some of the network elements (figure above) and are the following :
I) In the procambial position (central compartment), CK is initially available and there is an initial and sustained IAA input or selfupregulation.
This condition is supported by several lines of evidence. Also HB8, a marker of early vascular development that has been
found in preprocambial cells, is assumed to be initially present at this position. These conditions are not fixed, however. After the initial
configuration, all the members of the CK and IAA signaling pathways, as well as HB8, can change their states according to the logical rules.
II) In the xylem and phloem positions, it is assumed that no element is initially active except for the CK signaling pathway and TDIF, both in
the phloem position.The level of expression for a given node is represented by a discrete variable g and its value at a time t+1 depends on
the state of other components of the network (g1, g2, ..., gN) at a previous time unit. The state of every gene g therefore changes
according to:
gn(t+1)=Fn(gn1(t),gn2(t),..., gnk(t)) (1).
In this equation, gn1, gn2,…, gnk are the regulators of gene gn and Fn is a discrete function known as a logical rule (logical rules are
grounded in available experimental data, for example see slide 20). Given the logical rules, it is possible to follow the dynamics of the
network for any given initial configuration of the nodes expression state. One of the most important traits of dynamic models is the
existence of steady states in which the entire network enters into a selfsustained configuration of the nodes state. It is thought that in
developmental systems such self-sustained states correspond to particular cell types.
According to experimental evidence for the system under study, the hormone IAA, the peptide TDIF, and the microRNA MIR165/6 are able
to move among the cells. In the case of TDIF and MIR165/6, the mobility is defined as diffusion and is given by the following equation:
g(t+1)T[i]= H(g(t)[i]+ D (g(t)[i+1]+g(t)[i-1] – N(g(t)[i]))-b) (2),
where g(t)T[i] is the total amount of TDIF or MIR165 in cell (i). D is a parameter that determines the proportion of g that can move from any
cell to neighboring ones and is correlated to the diffusion rate of g. b is a constant corresponding to a degradation term. H is a step function
that converts the continuous values of g into a discrete variable that may attain values of 0, 1 or 2. N stands for the number of neighbors in
each cell. Boundary conditions are zero-flux. In the case of IAA, the mobility is defined as active transport dependent on the radial
localization of the PIN efflux transporters and is defined by the equation:
iaa(t+1)T[i]=Hiaa(iaa(t)[i]+Diaa(pin(t)[i+1])(iaa(t)[i+1])+Diaa(pin(t)[i-1])(iaa(t)[i-1])−N(Diaa)(pin(t)[i])(iaa(t)[i])−biaa) (3),
where Diaa is a parameter that determines the proportion of IAA that can be transported among cells. The transport depends on the
presence of IAA and PIN in the cells and biaa corresponds to a degradation term. As in equation 2, H is a step function that converts the
continuous values to discrete ones and N stands for the number of neighbors in each cell. Boundary conditions for IAA motion are also
zero-flux.
Using the logical rules, equations 1–3, and a broad range of parameter values (not shown here), it is possible fully to reproduce the results
and analyses reported in the following sections (see the figure above for the simulation time course).
35
36
Molecular Regulatory Networks Modeling
Benitez and Hejatko, submitted
□ The good model should be able to simulate reality
Another representation of the distinct expression profiles in the individual vascular
bundle compartments (phloem, procambium and xylem).
36
37
Molecular Regulatory Networks Modeling
Benitez and Hejatko, submitted
□ Simulation of mutants
37
38
38
 Definition of Systems Biology
 Tools
 Gene Ontology analysis
 Bayesian Networks
 Molecular/Gene Regulatory Networks Modeling
 Inferring Gene Regulatory Networks from Large Omics
Datasets
Outline
39
39
Systems Biology in
Cancer Research
40
40
miRNA/mRNA Profiling
Guo et al., Mol Med Reports, 2017
Hypertrophic scars (HS) area fibroproliferative disorder of the skin, which causes aesthetic and
functional impairment. However, the molecular pathogenesis of this disease remains largely
unknown and currently no efficient treatment exists. MicroRNAs (miRNAs) are involved in a variety
of pathophysiological processes, however the role of miRNAs in HS development remains unclear.
To investigate the miRNA expression signature of HS, microarray analysis was performed and 152
miRNAs were observed to be differentially expressed in HS tissue compared with normal skin
tissues. Of the miRNAs identified, miRNA-21 (miR-21) was significantly increased in HS tissues
and hypertrophic scar fibroblasts (HSFBs) as determined by reverse transcription-quantitative
polymerase chain reaction analysis. It was also observed that, when miR-21 in HSFBs was blocked
through use of an antagomir, the phenotype of fibrotic fibroblasts in vitro was reversed, as
demonstrated by growth inhibition, induction of apoptosis and suppressed expression of
fibrosis-associated genes collagen type I α 1 chain (COL1A1), COL1A2 and fibronectin.
Furthermore, miR-21 antagomir administration significantly reduced the severity of HS formation
and decreased collagen deposition in a rabbit ear HS model. The total scar area and scar elevation
index were calculated and were demonstrated to be significantly decreased in the treatment group
compared with control rabbits. These results indicated that the miR-21 antagomir has a therapeutic
effect on HS and suggests that targeting miRNAs may be a successful and novel therapeutic
strategy in the treatment of fibrotic diseases that are difficult to treat with existing methods.
miRNA expression signature profiling in hypertrophic scars (HS). (A) Volcano plot presenting
differentially expressed miRNAs between HS and paired (non-scar, obtained from donor sites
during scar resection) NS tissue. miRNA microarray expression profiling from three paired HS and
NS tissues was performed. Differentially expressed miRNAs were identified by fold change and a
P-value calculated using Student's t-test. The threshold set to identify up and downregulated genes
was a fold change ≥2 and P<0.05. Red dots indicate points-of-interest that exhibit large-magnitude
fold-changes (x-axis; log2 of the fold change) and high statistical significance (y-axis; -log10 of the
P-value). (B) Hierarchical clustering showing differentially expressed miRNAs from HS samples
compared with paired NS tissues. Each row represents one miRNA and each column represents
one tissue sample. The relative miRNA expression is depicted according to the color scale. Red
indicates upregulation and green indicates downregulation. N1-3 represents NS tissue samples,
whereas H1-3 represents HS tissue samples. The differentially expressed miRNAs were clearly
separated into clusters. miRNA, microRNA; hsa-miR, human microRNA; HS, hypertrophic scar;
NS, normal skin.
41
41 Benkova and Hejatko,Plant Mol Biol (2008)
proximalroot
meristem
distalroot
meristem
Inferring Gene
Regulatory Networks
42
Klidové centrum
Quiescent centre
Kolumela
Columella cell files
Iniciály kolumely
Columella initials
Epiderims
Epidermis
Kortex
Cortex
Endodermis
Endodermis
Iniciály stéle
Stele initials
proximalroot
meristem
distalroot
meristem
Postranní kořenová čepička
Lateral root cap
Iniciály epidermis
Epidermis initials
Iniciála endodermis a kortexu
Endodermis and cortex initial
Gene Regulatory
Networks
In the root, several functional and anatomical units could be recognized.
Along the longitudinal axis, the root meristem forms a distal root tip, including
stem cell niche, columella and lateral root cap, proximal meristem with a
population of rapidly dividing cells and elongation zone where cells leaving the
root meristem undergo rapid elongation and mature.
42
43 Birnbaum et al., Science, 2003
Gene Regulatory
Networks - GENIST
de Luis Balaguer et al., PNAS, 2017
 Inferring GRNs via GENIST
 GEne regulatory Network
Inference from SpatioTemporal
data algorithm
 Combining spatial- and timespecific
gene expression
profiles
43
44
Combining Large Omics
Datasets GENES
TISSUE/TIME
44
45
Gene Regulatory
Networks - GENIST
 Inferring GRNs via GENIST
 Clustering of genes
 Expression similarity under
various conditions/genetic
backgrounds, time points, …
 Inferring intra-cluster
connections
 Selection of potential
regulators and co-
regulators
 Based on the time
correlation in the change
of expression and/or user
specification
 Dynamic Bayesian
Network modeling
Haeseleer, Computational Biology, 2005
GENIST algorithm
The MATLAB source code for GENIST is publically available at
https://github.com/madeluis/GENIST.
For the detailed description of the procedure, see de Luis Balaguer et al., 2017,
SI (https://www.pnas.org/content/114/36/E7632/tab-figures-data)
45
46
 Inferring GRNs via GENIST
 Clustering of genes
 Expression similarity under
various conditions/genetic
backgrounds, time points, …
 Inferring intra-cluster
connections
 Selection of potential
regulators and co-
regulators
 Based on the time
correlation in the change
of expression and/or user
specification
 Dynamic Bayesian
Network modeling
de Luis Balaguer et al., PNAS, 2017
Gene Regulatory
Networks - GENIST
GENIST algorithm
The MATLAB source code for GENIST is publically available at
https://github.com/madeluis/GENIST.
For the detailed description of the procedure, see de Luis Balaguer et al., 2017,
SI (https://www.pnas.org/content/114/36/E7632/tab-figures-data)
GENIST block diagram. GENIST is implemented in MATLAB, and is composed of 
two consecutive steps, clustering and GRN inference. Clustering is performed 
based on a spatial dataset. Each resulting cluster is
independently processed by the GRN inference step, based on a temporal 
dataset.
46
47
Gene Regulatory
Networks - GENIST
de Luis Balaguer et al., PNAS, 2017
47
48 de Luis Balaguer et al., PNAS, 2017
Gene Regulatory
Networks - GENIST
48
49
Indirect PAN targets
feed-back loops
Gene Regulatory
Networks - GENISTMODEL PREDICTION
EXPERIMENTAL
VERIFICATION
PAN subnetwork in the QC inferred with the 12 developmental time points of the 
Arabidopsis root. (A) Optimal configuration (combination of signs— activation or 
repression—of the regulations that were inferred with undefined signs, which best fits 
the data in the simulations of the equations) of the subnetwork of PAN and its 
downstream targets. (B and C) Resulting expression values of PAN and its downstream 
targets, over time, after simulating the optimal configuration of the model. Simulations 
were run for 5 d and plots are shown until all factors reached steady states in the WT 
and pan mutant simulations. (B)Model simulated with the fitted equation parameters. 
(C)Model simulated with the PAN‐associated parameters set to zero to simulate a pan 
mutant situation. (D) Normalized expression values of PAN and its predicted 
downstream targets in Col‐0 wild type and in pan mutant. Statistically significant 
changes of expression between the mutant and the wild type, *q < 0.05.
In the WT simulation, all targets reached steady states by day 1 with subtle changes of 
expression during the transients (time length until expression values reach their steady 
states). On the contrary, the pan mutant simulation showed that EIN3 and WIP4 
presented high expression values during the transients and reached steady states at later 
stages (days 3 and 4, respectively). These delayed responses and initial activations of 
EIN3 and WIP4 reflect the prediction that these genes are indirectly affected by PAN. 
Further, the dynamics of our simulations depict that BRAVO, NTT, and WIP4 are, in our 
equations, connected through feedback loops. During the transient phase of the mutant 
simulation, NTT and BRAVO show an exponential decay, which is consistent with the 
prediction that they activate each other in the absence of PAN. However, their steady 
states are not immediately reached since they are activated by WIP4 and EIN3. 
Conversely, WIP4, which is repressed by a decaying NTT, shows high levels of expression.
49
50
50
 Definition of Systems Biology
 Tools
 Gene Ontology analysis
 Bayesian Networks
 Molecular/Gene Regulatory Networks Modeling
 Inferring Gene Regulatory Networks from Large Omics
Datasets
Summary
51
51
Discussion