Received: 5 September 2018
DOI: 10.1111/gtc.12655|

Revised: 4 November 2018|

Accepted: 5 November 2018

Genes to Cells

ORIGINAL ARTICLE

Reconstruction of global regulatory network from signaling tocellular functions using phosphoproteomic data
Kentaro Kawata1*

Yoko Tomizawa1||

Masaki Matsumoto8

Katsuyuki Yugi1,2,3,4*

Masashi Fujii1,6

|

|

|

Atsushi Hatano1

Shinsuke Uda7

Keiichi I. Nakayama8

|

|

|

Toshiya Kokaji5

Hiroyuki Kubota7

Shinya Kuroda1,5,9

|

|

1

Department of Biological Sciences, Graduate School of Science, University of Tokyo, Bunkyo‐ku, Japan

2

YCI Laboratory for Trans‐Omics, Young Chief Investigator Program, RIKEN Center for Integrative Medical Science, Yokohama, Japan

3

Institute for Advanced Biosciences, Keio University, Fujisawa, Japan

4

PRESTO, Japan Science and Technology Agency, Yokohama, Japan

5

Department of Computational Biology and Medical Sciences, Graduate School of Frontier Sciences, University of Tokyo, Kashiwa, Japan

6

Molecular Genetics Research Laboratory, Graduate School of Science, University of Tokyo, Bunkyo‐ku, Japan

7

Division of Integrated Omics, Research Center for Transomics Medicine, Medical Institute of Bioregulation, Kyushu University, Fukuoka, Japan

8

Department of Molecular and Cellular Biology, Medical Institute of Bioregulation, Kyushu University, Fukuoka, Japan

9

Core Research for Evolutional Science and Technology (CREST), Japan Science and Technology Agency, Bunkyo‐ku, Japan

Correspondence
Shinya Kuroda, Department of Biological
Sciences, Graduate School of Science,
University of Tokyo, Bunkyo‐ku, Japan.
Email: skuroda@bs.s.u-tokyo.ac.jp
Present Address
Kentaro Kawata, Isotope Science
Center, University of Tokyo, Bunkyo‐ku,
Japan
Funding information
Japan Science and Technology Corporation,
Grant/Award Number: JPMJPR1538;
Grant‐in‐Aid for Scientific Research on
Innovative Areas, Grant/Award Number:
18H04801; Core Research for Evolutional
Science and Technology, Grant/Award
Number: JPMJCR12W3; Japan Society for
the Promotion of Science, Grant/Award
Number: 16K12508, 16K12508, 16H06577,
17H06299, 17H06300, 17H06301,
17H06306, 18H03979, JP15H05582
and JP18H05431; National Institute of
Genetics (NIG); Research Organization ofInformation and Systems (ROIS)
Communicated by: Kozo Kaibuchi

Abstract
Cellular signaling regulates various cellular functions via protein phosphorylation.
Phosphoproteomic data potentially include information for a global regulatory network from signaling to cellular functions, but a procedure to reconstruct this network
using such data has yet to be established. In this paper, we provide a procedure toreconstruct a global regulatory network from signaling to cellular functions from
phosphoproteomic data by integrating prior knowledge of cellular functions and inference of the kinase–substrate relationships (KSRs). We used phosphoproteomic
data from insulin‐stimulated Fao hepatoma cells and identified protein phosphorylation regulated by insulin specifically over‐represented in cellular functions in the
KEGG database. We inferred kinases for protein phosphorylation by KSRs, and connected the kinases in the insulin signaling layer to the phosphorylated proteins in thecellular functions, revealing that the insulin signal is selectively transmitted via thePi3k‐Akt and Erk signaling pathways to cellular adhesions and RNA maturation,
respectively. Thus, we provide a method to reconstruct global regulatory networkfrom signaling to cellular functions based on phosphoproteomic data.
KEYWORDS
cellular functions, insulin signaling, network integration, phosphoproteomic data, phosphorylation,systems biology, Trans‐omics*

These authors contributed equally to this work.

Genes Cells. 2018;1–12.

wileyonlinelibrary.com/journal/gtc© 2018 Molecular Biology Society of Japan and John
Wiley & Sons Australia, Ltd

| 1

2

|   

1

|

Genes to Cells
IN TRO D U C T ION

Protein phosphorylation plays a pivotal role in signaling‐dependent regulation of various cellular functions. The recent
development of proteomic technologies has enabled us toquantify abundance of protein phosphorylation comprehensively (Humphrey, James, & Mann, 2015; Olsen et al., 2006).
Phosphoproteomic data potentially include information about
global regulatory networks from signaling to cellular functions via protein phosphorylation. Prior knowledge of cellularfunctions in databases such as KEGG which is a manually curated database for biological pathway and functions of genes(Kanehisa, Furumichi, Tanabe, Sato, & Morishima, 2017;
Kanehisa, Goto, Sato, Furumichi, & Tanabe, 2012), reactome
which is a manually curated database for biological pathways
and processes of signaling and metabolic molecules (Croft etal., 2014), and Gene Ontology (GO) which is a set of functionalterms related to biological process and localization of proteinsrepresented as directed acyclic graph (Ashburner et al., 2000)
has been used to identify specific cellular functions of protein
phosphorylation from phosphoproteomic data. Several tools
inferring kinase–substrate relationship (KSR) based on similarity of recognition motifs have also been developed, suchas NetPhorest (Miller et al., 2008), KinomeXplorer (Horn etal., 2014), KinasePhos (Huang, Lee, Tzeng, & Horng, 2005)
and iGPS (Song et al., 2012), and using such methods, association of kinases with protein phosphorylation can be inferredfrom phosphoproteomic data (Emdal et al., 2015; Giansantiet al., 2015; Kanshin, Kubiniok, Thattikota, D’Amours, &Thibault, 2015; Patella et al., 2015; Sacco et al., 2016; Yugi
et al., 2014). Biological interpretation of phosphoproteomic
data based on prior knowledge has addressed specific cellularfunctions regulated by protein phosphorylation, but not thesignaling‐dependent regulation of such cellular functions.
KSR addresses the regulation of substrates by kinases in thesignaling pathway, but not regulation of the cellular functions.
Integrating these two approaches may enable the reconstruction of the global regulatory network from signaling to cellular functions based on phosphoproteomic data, but a method
for this has not thus far been established.
Insulin regulates diverse cellular functions throughphosphorylation of substrate proteins by protein kinases in
the insulin signaling pathway, including Akt and extracellular signal‐regulated kinase (Erk) (Lizcano & Alessi, 2002;
Saltiel & Kahn, 2001); insulin‐regulated cellular functions
include cell adhesion, cytoskeletal organization, RNA transport and RNA splicing (Hartmann et al., 2009; Reiss et al.,
2001; Tsakiridis et al., 1999; Wolf et al., 2013). It has beenreported that mutation of substrate proteins for Akt and inhibition of phosphoinositide 3‐kinase (Pi3k) reduce actin
reorganization and intracellular adhesion in response to insulin in cancer cells (Reiss et al., 2001; Tsakiridis et al.,
1999; Wolf et al., 2013). Erk, a member of insulin signaling,

KAWATA et al.

has also been reported to regulate multiple alternative splicing factors, as well as playing a role in transcriptional regulation (Al‐Ayoubi, Zheng, Liu, Bai, & Eblen, 2012; Boulton
et al., 1991; Eblen, 2018; Matter, Herrlich, & König, 2002).
Thus, insulin regulates various cellular functions by phosphorylating multiple substrate proteins via various protein
kinases. However, the global regulatory network of insulin
action via protein phosphorylation has not thus far been
elucidated (Buescher et al., 2012; Chiappino‐Pepe, Pandey,
Ataman, & Hatzimanikatis, 2017; Hatzimanikatis & Saez‐
Rodriguez, 2015; Hyduke, Lewis, & Palsson, 2013; Joyce &Palsson, 2006; Palsson & Zengler, 2010). Phosphoproteomic
studies of insulin signaling have been performed with
Drosophila cells (Friedman et al., 2011; Vinayagam et al.,2016), 3T3‐L1 mouse adipocytes (Humphrey et al., 2013),
mouse hepatoma cells (Monetti, Nagaraj, Sharma, & Mann,
2011), mouse brown preadipocytes (Krüger et al., 2008),
rat primary hepatocytes (Zhang, Zhang, & Yu, 2017), rat
Fao hepatoma cells (Yugi et al., 2014) and mouse hepatoma (Humphrey, Azimifar, & Mann, 2015). Among these
phosphoproteomic studies, the cellular functions of insulin‐dependent protein phosphorylation have been estimatedusing databases of prior knowledge such as GO (Friedman
et al., 2011; Humphrey, Azimifar, et al., 2015; Vinayagam
et al., 2016; Zhang et al., 2017), whereas insulin‐dependent
regulatory networks of protein phosphorylation by kinases
have been estimated using KSR inference (Humphrey et al.,2013; Miller et al., 2008). However, none of these studies
integrated both prior knowledge of cellular functions fromthe databases and KSR inference, so the global regulatory
network from insulin signaling to cellular functions has notthus far been reconstructed.
This study integrates prior knowledge of cellular functions from the above‐mentioned databases and KSR inference to provide a procedure to reconstruct the global
regulatory network from signal to cellular function basedon phosphoproteomic data. First, we identified the cellular functions of protein phosphorylation from phosphoproteomic data by over‐representation analysis using the
KEGG database. Next, we inferred the protein kinases
regulating protein phosphorylation in the cellular functions
by KSR using NetPhorest, along with the selective use ofkinases for each cellular function. Finally, we mapped theinferred protein kinases in the insulin signaling pathway
in the KEGG database, by connecting the protein kinasesin the signaling pathway to the protein phosphorylationin the cellular functions layer and reconstructing global
regulatory network from insulin signaling to the cellularfunctions via protein phosphorylation. The global regulatory network revealed that the insulin signal is selectivelytransmitted via the Pi3k‐Akt pathway and the Erk pathway
to cellular adhesions and RNA maturation such as RNA
transport or splicing, respectively.

KAWATA et al.

(a)

(b)

79 (0.93%)
1,421(16.74%)

6,989 (82.33%)

2

|

pSer
pThr
pTyr

515 (6.50%)

26 (0.33%)
1 (0.01%)| 3

Singly phosphorylated
Doubly phosphorylated
Triply phosphorylated
Quadruply phosphorylated

7,377 (93.16%)

250

Proteins

200

IRpPs

150
100
50
0

tj
un
ct
Sp
io
n
lic
eo
Lo
s
o
RN
ng
m
e
ev
A
tra
ity
Lo
ns
re
ng
po
gu
ev
rt
la
ity
tin
- m reggpul ulaat
tip ti
h.
le ng
sp p
Fo
a
e
ca cie th.
la s
Ad
dh
he
es
re
io
n
ns
ju
nc
t
io
Ph
n
A
o
si sph pop
gn atoal tidin yli sis
g no
sy s
st ito
e l
G
Re
ap m
gu
ju
n
la
c
cy tion tion
to o
sk f ael c
et tinon

Numbers of proteins

(c)

Ti
gh

F I G U R E 1 Characterization ofphosphoproteomic data. Distribution ofphosphorylation sites by (a) amino acidresidues and (b) number of phosphorylatedsites per peptide. Numbers and parentheses
indicate the number of phosphorylationsites or phosphopeptides, and percentage
relative to the total phosphorylation sites orphosphopeptides, respectively. (c) Numberof proteins in the KEGG cellular functional
pathways. The light gray bars indicate thetotal number of proteins, and the dark graybars indicate the number of IRpP containedin each signaling pathway  Genes to Cells

R E S U LTS

2.1 | Characterization of thephosphoproteomic data
Previously, we obtained phosphoproteomic data of rat Fao
hepatoma cells in acute insulin action (<60 min; Yugi et al.,2014). We used the 199 phosphopeptides from the 49 metabolic enzymes in our previous work; the remaining phosphopeptides have yet to be analyzed. For this analysis, we use all
7,929 phosphopeptides from 3,468 proteins; the phosphopeptides included 6,989 phosphoserine (pSer), 1,421 phosphothreonine (pThr) and 79 phosphotyrosine (pTyr) (Figure 1a).
Singly phosphorylated peptides represented 93% of the 7,929
phosphopeptides (Figure 1b). The distribution of phosphorylation sites by amino acid was consistent with the previous
studies (Humphrey et al., 2013; Humphrey, Azimifar, et al.,2015; Olsen et al., 2006, 2010).

2.2 | Identification of cellular functions
regulated by insulin
To identify the cellular functions regulated by insulin via protein phosphorylation, we identified quantitatively changed
phosphopeptides from the phosphoproteomic data. We identified 3,288 quantitatively changed phosphopeptides exhibiting a change in phosphorylation intensity greater than
a 1.5‐fold increase or less than a 0.67‐fold decrease at one
or more time points in response to insulin stimulation; thesesites were present on 1,947 proteins. Hereafter, we define a

protein with at least one quantitatively changed phosphopeptide as an insulin‐responsive phosphoprotein (IRpP).
To determine cellular functions related to the 1,947 IRpPs,
we performed pathway over‐representation analysis using
KEGG pathways (that is, those pathways included in the
KEGG database). We extracted pathways related to various
cellular functions from KEGG (hereafter cellular functional
pathways), except for those related to signaling and those thatfunction specifically in tissues other than liver (see Section
4). The analysis with KEGG pathways related to signaling
was performed in our previous study (Kawata et al., 2018).
The pathways that function specifically in tissues other thanliver were excluded because the phosphoproteomic data arefrom rat Fao hepatoma cells, which are a liver cell line. IRpPs
were significantly over‐represented in 11 out of 128 cellularfunctional pathways, including those regulating cell adhesion(Figure 1c), such as Tight junction (rno04530) and Gap junction (rno04540), and those regulating RNA maturation, such
as RNA transport (rno03013) and Spliceosome (rno03040).
We extracted the 11 pathways in which IRpPs were significantly over‐represented and defined them as the cellularfunctional layer (Table 1, Supporting Information Table S1).

2.3 | Inference of protein kinases forprotein phosphorylation in the cellularfunctional pathways
To infer the specific protein kinases regulating the IRpPs
in the cellular functional pathways, we inferred protein kinases recognizing quantitatively changed phosphopeptides as

4

|   

KAWATA et al.

Genes to Cells

TABLE 1

Cellular functional pathways in which quantitatively
changed phosphopeptides were over‐represented
Pathway nameOdds ratio

p value

FDR−9

Tight junction

4.11

2.12 × 10

2.41 × 10−7

Spliceosome

3.40

7.64 × 10−7

2.85 × 10−5

Longevity regulating pathway

3.24

−5

7.61 × 10

1.16 × 10−3

RNA transport

2.59

7.17 × 10−5

1.16 × 10−3

Longevity regulating pathway—multiple species

3.73

1.30 × 10−4

1.77 × 10−3

Adherens junction

3.16

5.45 × 10−4

4.92 × 10−3

Focal adhesion

2.17

−4

5.07 × 10

4.92 × 10−3

Apoptosis

2.11

4.43 × 10−3

3.15 × 10−2

Phosphatidylinositol
signaling system

2.30

6.94 × 10−3

4.72 × 10−2

Gap junction

2.33

8.20 × 10−3

5.33 × 10−2

Regulation of actin
cytoskeleton

1.74

1.19 × 10−2

7.40 × 10−2

substrates in the cellular functional pathways, using the KSR
prediction tool NetPhorest (Horn et al., 2014; Miller et al.,2008). The kinase classifier with the largest posterior probability value for each phosphopeptide was considered the
responsible protein kinase. We inferred kinase–substrate relationships between a total of 26 responsible protein kinases
and 489 out of 492 quantitatively changed phosphopeptides
from the 216 IRpPs in the 11 cellular functional pathways(Supporting Information Table S2).

2.4 | Selective use of protein kinases
according to cellular functions
To examine the selective use of the protein kinases regulating specific cellular functions, we calculated the occurrence rate of the responsible kinases of the phosphopeptidesof IRpPs in each cellular functional pathway (Supporting
Information Table S3). We performed cluster analysis using
the occurrence rates, which resulted in two clusters (Figure2a): Cluster 1 included pathways regulating cell adhesion

and the cytoskeleton, such as Tight junction (rno04530),
Gap junction (rno04540) and Regulation of actin cytoskeleton (rno4810) (Figure 2a, blue); Cluster 2 included pathways regulating RNA maturation, such as RNA transport(rno03013) and Spliceosome (rno03040) (Figure 2a, red).
To investigate the protein kinase regulating each cellular function, we compared the motif logos of amino acidsequences (Colaert, Helsens, Martens, Vandekerckhove, &Gevaert, 2009) of the quantitatively changed phosphopeptides of IRpPs in each cellular functional pathway (Figure2b). For many of the pathways in Cluster 1, arginine residue(R) at the −3 position was emphasized, whereas proline residue (P) at the +1 position was emphasized for many pathways
in Cluster 2. R at the −3 position and P at the +1 position
were known as target sequences of basophilic kinases suchas Akt and proline‐directed kinases such as Erk, respectively.
The results of cluster analysis and motif logos indicate thatcell adhesion and the cytoskeleton were selectively regulated
mainly by Akt, whereas RNA maturation was selectively regulated mainly by Erk.
To confirm these findings, the representation of theamino acid residues at each position of the phosphopeptides(Workman et al., 2005) included in each cluster was tested
statistically (Figure 2c). In phosphopeptides of IRpPs on thepathways in the cellular functional layer, P at the +1 position
and R at the −3 position were significantly over‐represented(p < 0.05) in comparison with the amino acid compositionof rat proteins, indicating insulin regulates cellular functions
via basophilic kinases such as Akt and proline‐directed kinases such as Erk. In the phosphopeptides in Cluster 1, only
a valine residue (V) at the +1 position was over‐represented
and P at the +1 position was significantly under‐represented
in comparison with all phosphopeptides in the phosphoproteomic data; contrary to the expectations, R at the −3 position
was not significantly over‐represented. For the phosphopeptides in Cluster 2, P at the +1 position was significantly
over‐represented. These results suggest selective regulationof cellular functions by two distinct type of protein kinases:
Akt regulates both of the pathways in Cluster 1 (such as celladhesion and the cytoskeleton) and those in Cluster 2 (suchas RNA maturation), whereas Erk regulates the pathways in
Cluster 2 specifically.

F I G U R E 2 Protein kinase selectivity for protein phosphorylation in the cellular functional pathways. (a) Clustering of cellular functional
pathways by occurrence rates for each kinase class. Green indicates the larger occurrence rate of the kinase class. (b) Motif logos of amino acidsequences of phosphopeptides in each of the indicated pathways. The height of each letter at each position is scaled relative to the information
content, reflecting the frequency of the corresponding amino acid. Blue pathway names correspond to those in Cluster 1; red names correspond to
those in Cluster 2. The number in parentheses attached to each pathway name represents the number of quantitatively changed peptides includedin the pathway. Note that statistical differences of frequencies from any reference sequences are not considered for these motif logos. (c) Statisticalamino acid motifs for phosphopeptides of IRpPs included in the cellular functional layer (left), Cluster 1 (middle) and Cluster 2 (right) are shownas logo plots. The characters above and below each horizontal line indicate amino acids showing significant over‐representation and under‐
representation (p < 0.05), respectively. The height of the letter representing an amino acid at each position reflects the difference in the frequency of
its occurrence in the sets of phosphopeptides and reference proteins

KAWATA et al.

  Genes to Cells
Kinase class(a)

Atypical

STE

CMGC

AGC

TK

| 5

Occurrence rate ofkinase class

CAMK

0.0 0.05 0.1 0.15 0.2 0.25

Cluster 1

Apoptosis

Cellular functional pathways

Adherens junction
Tight junction
Phosphatidylinositol signaling system
Focal adhesion
Regulation of actin cytoskeleton

Cluster 2

Longevity regulating pathway - multiple species
Longevity regulating pathway
Gap junction
Spliceosome
RNA transport

0.4

0.2

0.0

AC
T
A T PD
FL R2_
T3 AC
M HK
_A _g
_C T
TR ro
SF R2
u
1R B_
C _g pT_K G K1 rou
it_ Fb _g p
PD R2 rou
G _g pF r
EG R_ oup
FR grou
SL _gr p
K_ oup
C
Y
g
D
M
SK ro
AP K2
u
P _g p
K3 _C
M AK rou_M DK
AP _ g p
3
AP _2K rouK1 CD__M K1 CK gro p_
AP CD 2_ up
K7 K5 gro
_ N _ g up
LK rou_p3 gro p
G 8_ up
SK gr
3 ou
RC _g pK_ roup
CL gr
K ou
PK _gr p
DM A_ oup
PK groG _ g upRK rou
PK _gr p
B ou
PK _gr p
DA C_ oup
PK gro_ g uprou
PP p
Ca KD im
M _gr 2
KI ou
I_
gr p
ou
p

0.6

Adherens junction (56)

Apoptosis (49)

bits

(b)
3.0

3.0

1.5

1.5

0.0

N -5 -4 -3 -2 -1 0 +1 +2 +3 +4 +5 C

0.0

bits

3.0

1.5

0.0

N -5 -4 -3 -2 -1 0 +1 +2 +3 +4 +5 C

bits

1.0
0.0

N -5 -4 -3 -2 -1 0 +1 +2 +3 +4 +5 C

Relative location
(c)

2.0

1.5

N -5 -4 -3 -2 -1 0 +1 +2 +3 +4 +5 C

0.0

signaling system (29)1.0
N -5 -4 -3 -2 -1 0 +1 +2 +3 +4 +5 C

0.0

3.0

3.0

1.5

1.5

1.5

0.0

N -5 -4 -3 -2 -1 0 +1 +2 +3 +4 +5 C

0.0

RNA transport (94)

Spliceosome (47)
3.0

3.0

1.5

1.5

0.0

N -5 -4 -3 -2 -1 0 +1 +2 +3 +4 +5 C

Relative location

Cellular functional pathways

N -5 -4 -3 -2 -1 0 +1 +2 +3 +4 +5 C

0.0

N -5 -4 -3 -2 -1 0 +1 +2 +3 +4 +5 C

Longevity regulating
pathway (54)Longevity regulating pathway- multiple species (40)3.0

Gap junction (40)
2.0

3.0

Regulation of actin
cytoskeleton (80)

Focal adhesion (74)

Phosphatidylinositol

Tight junction (130)0.0

N -5 -4 -3 -2 -1 0 +1 +2 +3 +4 +5 C

Relative location

N -5 -4 -3 -2 -1 0 +1 +2 +3 +4 +5 C

Relative location

Cluster 1

Cluster 2

11

25

52

5

12

0 -5 -4 -3 -2 -1 0 1 2 3 4 5

0 -5 -4 -3 -2 -1 0 1 2 3 4 5

0 -5 -4 -3 -2 -1 0 1 2 3 4 5

p value = 0.05
% difference

105

-52

-5

-12

-105

-11

-25

6

|   

KAWATA et al.

Genes to Cells

2.5 | Reconstruction of global regulatory
network from signaling to the cellular functions
Finally, we reconstructed the global regulatory networkfrom signaling to the cellular functional layer via protein
phosphorylation by connecting the inferred responsible
kinases in the insulin signaling pathway to IRpPs in thecellular functional layer. As the insulin signaling pathway, we used insulin signaling pathway (rno04910) in the
KEGG database (Figure 3a). Among the 61 molecules on
the signaling pathway, 12 are registered as protein kinases
in NetPhorest. Among the 26 kinase groups inferred in
Figure 2a, Akt in the PKB_group, Erk1/2 in the MAPK3_
MAPK1_MAPK7_NLK_group, Pka in the PKA_group,Pkc in the PKC_group, Mek1/2 in the MAP2K_group andGsk3β in the GSK3_group were included in the insulin
signaling pathway (Figure 3a, red). Mek1/2 was the only
kinase in MAP2K_group involved in insulin signaling pathway. Given that Mek1/2 is a specific kinase for Erk family(Pearson et al., 2001), other kinases in the MAP2K_group
not involved in insulin signaling pathway, rather than
Mek1/2, may be the responsible kinase for the IRpPs regulated by MAP2K_group.
We connected the protein kinases in the insulin signaling pathway to the IRpPs in the cellular functional layer, resulting in a reconstructed global regulatory network (Figure3b). Among the 489 phosphopeptides of IRpPs in the cellularfunctional layer for which responsible kinases were inferred,
193 phosphopeptides of IRpPs were connected to the signaling pathway via the kinases (Figure 3b).
To examine the functional selectivity of the protein kinases, we counted the number of phosphopeptides of IRpPs
in each cluster connected to the protein kinases (Table 2). Of
the 193 phosphopeptides of IRpPs in the cellular functional
layer, 122 peptides (63.2%) and 93 peptides (48.2%) were
included in Clusters 1 and 2, respectively; 22 phosphopeptides (11.4%) were included in both clusters. Kinases in thePi3k‐Akt signaling pathway such as Pka, Pkc and Akt wereconnected to more phosphopeptides in Cluster 1 than Cluster
2, whereas the kinases in the Erk signaling pathway—such as
Mek1/2 and Erk1/2—were connected to more phosphopeptides in Cluster 2 (Figure 3c). This indicates that the Pik3‐
Akt and Erk signaling pathways selectively regulate different
cellular functions via protein phosphorylation; the former
mainly regulates cell adhesion such as the adherens junction,and the latter mainly regulates RNA maturation such as RNA
splicing and RNA transport.

3|

DISCUSSION

We analyzed the phosphoproteomic data of insulin‐stimulated Fao hepatoma cells by integrating prior knowledge ofcellular functions and KSR inference and reconstructed a
global regulatory network from signaling to cellular function. This reconstructed network of insulin action revealed
that insulin signals are selectively transmitted via the Pi3k‐
Akt and Erk signaling pathways to cellular adhesions andRNA maturation, respectively. This result is consistent withthe earlier finding that Akt and Erk play important roles ininsulin signaling pathway in adipocytes using phosphoproteomic data (Humphrey et al., 2013). Blood insulin showeda transient high concentration during the fed state (induced
insulin secretion), whereas a sustained low concentrationof blood insulin (basal insulin secretion) was maintained
during the fasting state (Lindsay et al., 2003; Polonsky,
Given, & Van Cauter, 1988). We have reported that highand low doses of insulin selectively regulate different signaling molecules in the liver (Kubota et al., 2012). We also
found that high and low doses of insulin regulate different
cell functions through Erk and Akt, in Fao rat hepatoma
cells and the rat liver (Kawata et al., 2018). Our finding
suggests that high and low doses of insulin selectively regulate cellular adhesions via the Pi3k‐Akt signaling pathway and RNA maturation via the Erk signaling pathway,
respectively, in vivo. Elucidating such selective regulationin vivo is a subject for future study.
Previous studies combining phosphoproteomic data and
NetPhorest indicated that Akt has an important role for cellular regulation in response to insulin (Alli Shaik et al., 2016;
Vinayagam et al., 2016), consistent with this study. These studies did not focus on the kinase selectivity of the regulated cellular functions, and this study provides the first evidence for the
kinase selectivity of cellular functions regulated by insulin. The
similarities in phosphorylation site consensus motifs shared bymany protein kinases may also have resulted in inaccurate predictions of KSR. Although approximately 500 protein kinases
have been reported in humans (Manning, Whyte, Martinez,
Hunter, & Sudarsanam, 2002), only approximately 200 kinases
are registered in NetPhorest, which is the KSR prediction tool

F I G U R E 3 Reconstruction of the global regulatory network from the insulin signaling pathway to the cellular functional layer. (a) Insulin
signaling pathway based on KEGG database used as the signaling pathway in the global regulatory network. Red and blue indicate kinases
connected to the phosphopeptides of the IRpPs, and other kinases, respectively. Pi3k‐Akt and Erk signaling pathways are indicated by gray frames.
(b) Constructed global regulatory network. The arrows (from the top to the bottom layer) indicate phosphorylation of the quantitatively changed
phosphopeptides by the responsible protein kinases. (c) Number of quantitative changed phosphopeptides inferred as substrates for the kinases on
the signaling pathway. Red and blue bars indicate the number of phosphopeptides on only the pathways in Clusters 1 and 2, respectively; the green
bar indicates the number of phosphopeptides on pathways in both clusters

KAWATA et al.

Genes to Cells(a)GlucoseFlotillin

Cap

CrkII

Cbl

Grf2

Aps

Exo70

Tc10

Glut4

Cip4/2

Ampk

Ptp1b

Srebp1

Pyk
G6pc

Foxo1

Pkc

Ship2

Acc
Fas

Pgc1a
Socs

  

Gk

Fbp
Pepck

Insulin

InsR

Irs

Pi3k

Pdk1/2

PIP3

Akt

Gsk3β

Tsc2

Pp1

Tsc1

Jnk
Ikk

mTor

Lar

Gys
Phk

Pde3b

cAMP

Pka

Glycogen

Pyg

Hsl

Rheb
Bad

Raptor

S6k

Pi3k-Akt signaling
pathway

S6
eIf4e

eIf4ebp1
Mnk

Shc

Grb2

Sos

Mek1/2

Raf

Ras

Erk1/2

Elk1

Erk signaling pathway(b)

Insulin signaling pathway

InsulinPkc
Ampk

GlucoseFlotillin

CapAps

6 responsible
protein kinases

CrkII

Cbl

Exo70

Tc10

Grf2

Glut4

Cip4/2

Srebp1

Pyk

Pgc1a
Socs
InsR

Irs

Ptp1b
Pi3k

Pkc

Ship2

Foxo1

Pp1

Tsc2
Tsc1

Mek1/2

mTor
S6k

Grb2

Sos

Ras

Raf

Mek1/2

Mnk

Erk1/2

Pka

Glycogen

Gys
Phk

cAMP

Gk

Pka

Pyg
Hsl

Erk1/2Bad

S6

eIf4ebp1
Shc

Pde3b

Rheb

Raptor

G6pc
Fbp
Pepck

Gsk3β

Akt

Pdk1/2

PIP3

Jnk
Ikk

Lar

Akt
Gs k3β

Acc
Fas

eIf4e
Elk1

a. Tight junction
b. Focal adhesionc. Regulation ofactin cytoskeleton
d. Apoptosis
e. Adherens junction
f. Phosphatidylinositol
signaling system

112 IRpPs
192 quantitatively
changed phosphopeptides

h
k
i

g

Longevity regulating pathway

RNA transport

Longevity regulating pathway- multiple speciesjSpliceosomef

Gap junctionedApoptosis

Adherens junction

Phosphatidylinositol
signaling system

c

b

a

Tight junction

Focal adhesion

k. RNA transport

Numbers of peptides

60
Cluster 1
Both Clusters
Cluster 2

50
40
30
20
10
0
a
Pk

(P

KA

(
1/2
Erk

MA

P

)

p)
p)
p)
p)
p)
rou
rou
rou
rou
rou
g
g
g
g
g
_
_
3_
K_
C_
2K
KB
SK
NL
PK
AP
t (P
(G
7_
c(
k
M
k
K
(
β
A
P
P
/2
k3
k1
MA
Gs
1_
Me
K
AP
_M
3
K

p
rou
_g

i. Longevity regulating pathway- multiple species
j. Spliceosome

Regulation ofactin cytoskeleton

Cellular functional layer(c)g. Gap junction
h. Longevity regulating pathway

| 7

8

|   

KAWATA et al.

Genes to Cells

TABLE 2

Number of quantitatively changed phosphopeptides
inferred as substrates for kinases on the signaling pathway
KinaseCluster 1

Cluster 2

Both clusters

Total

Pka

48

6

2

56

Pkc

28

15

7

50

Erk1/2

5

35

2

42

Akt

17

8

3

28

Mek1/2

1

5

8

14

Gsk3β

1

2

0

3

Total

100

71

22

193

with the largest number of registered protein kinases. Moreover,the phosphoproteomic data in this study obtained with shotgun
approach include missing KSRs and nonfunctional phosphorylation. We select the pathway in which IRpPs are concentratedusing statistical over‐representation analysis to reduce the influence of missing and nonfunctional phosphorylation on thelandscape of insulin action. Although KSR prediction can findnovel pathways between the signaling pathway and cellularfunction, it cannot find novel cellular functions. In this study,
we used KEGG pathways to provide prior knowledge of cellularfunctions because they clearly distinguish signaling and cellularfunctions. In addition to KEGG pathways, GO (Ashburner et al.,2000) and pathway information included in reactome (Croft etal., 2014) are used to provide established pathway information.
Depending on cellular function, the precision and completenessof prior knowledge may be quite different; it is possible thatover‐representation may not have resulted in the original cellularfunction. Despite these limitations, we obtained novel molecular insight into how cells interpret insulin stimulation throughthe integration of cellular functions available in the databases
mentioned above to the signaling identified by KSR prediction.
Since metabolism did not indicate significant over‐representation in the pathway analysis, it is not included in thetargeted cellular functions. However, metabolic regulationis one of the well‐known functions of insulin (Lizcano &Alessi, 2002; Saltiel & Kahn, 2001), and we have previously
reconstructed a global regulatory network involved in metabolism (Yugi et al., 2014). We revealed the large contributionof Akt and Erk in regulating cellular functions via phosphorylation in response to acute insulin action. We propose thatthese reconstruction methods be applied to explore cellularresponses to other stimuli.

4 | EX P ER IME N TA L
PRO C E D U R E S
4.1|

Phosphoproteomic data

In this study, we used published datasets of the quantitative
phosphoproteomic data of a time series of insulin stimulationof Fao cells (RRID:CVCL_0269, male) that we measured
previously (Yugi et al., 2014). Rat hepatoma Fao cells wereseeded at a density of 3 × 106 cells per dish on 6‐cm dishes(Corning) and cultured in RPMI 1640 supplemented with10% (v/v) fetal bovine serum at 37°C under 5% CO2 for
2 days before deprivation of serum (starvation). The cells
were washed twice with phosphate‐buffered saline (PBS) andstarved for 16 hr in serum‐free medium including 0.01 nM
insulin (Sigma‐Aldrich) and 10 nM dexamethasone (Wako),which increases the expression of gluconeogenesis genes suchas G6pase and Pck1 (Lange et al., 1994). We continuously
added 0.01 nM insulin before the stimulation, and 0.01 nM
insulin was present throughout the experiments unless otherwise specified to mimic in vivo basal secretion during fasting (Polonsky et al., 1988). The medium was changed at 4
and 2 hr before the stimulation. Cells were stimulated with
1 nM insulin and collected at 0, 2, 5, 10, 30, 45 and 60 min
after stimulation. Cell lysate digested with LysC and trypsin
was subjected to Fe‐IMAC and iTRAQ labeling for the enrichment of phosphopeptides and quantification by mass
spectrometry. All samples were analyzed with a QSTAR
Elite (AB Sciex) instrument equipped with a Paradigm MS4
HPLC pump and HTC‐PAL autosampler (CTC Analytics
AG). The peak lists were generated using Analyst Mascot.dll
v1.6b27 (AB SCIEX). A MASCOT search was performed
with the following parameter settings: trypsin as the enzyme
used; the allowed number of missed cleavages as 2; iTRAQ
label at the NH2‐terminus, Lys and carbamidomethylationof Cys as fixed modifications; oxidized Met, iTRAQ label
on Tyr, pyroglutamination of NH2‐terminal Glu or Gln, andphosphorylation on Ser, Thr and Tyr as variable modifications; precursor mass tolerance as 100 ppm; and tolerance of
MS/MS ions as 0.2 Da. Assigned rank 1 peptide sequences(MASCOT score >20) were extracted. Evaluation of phosphorylation sites was performed at a post‐MASCOT searchwith in‐house script.

4.2 | Identification of differential
phosphorylation
Because the phosphoproteomic data consisted of two different time series from two separate experiments (0, 5, 10
and 45 min and 2, 10, 30 and 60 min), some of the phosphopeptides were identified and quantified in data from only
one of the time series. Therefore, a fold change of phosphorylation intensity was calculated as a ratio of the phosphorylation intensity at each time point to the phosphorylation
intensity at t = 0 or 2 min. A phosphopeptide with a phosphorylation intensity greater than a 1.5‐fold increase or lessthan a 0.67‐fold decrease at one or more time points was defined as a quantitatively changed phosphopeptide. Proteins
including one or more quantitatively changed phosphopeptides were defined as insulin‐responsive phosphoproteins

KAWATA et al.

(IRpPs). The detailed procedure of identification of differential phosphorylation can be found in Yugi et al (Yugi etal., 2014).

4.3 | Construction of the cellularfunctional layer
We selected all rat pathways from the KEGG databasefor cellular functions with the following exceptions: (a)signaling pathways (43 pathways); (b) global pathways(rno01XXX) including Metabolic pathways (rno01100)(9 pathways); (c) disease‐related pathways (rno05XXX)(63 pathways); (d) pathways that work in specific tissues
other than liver or are not regulated by insulin signaling in
the liver (52 pathways); and (e) pathways that include thecharacter string of “diabetes,” “NAFLD” or “Insulin resistance”; in their names (8 pathways). The reasons why we
excluded the above five kinds of pathways are as follows.
The signaling pathways were excluded from cellular functional pathways because these were analyzed in the previous study (Kawata et al., 2018). Global and disease‐related
pathways were excluded because of their redundancy withother cellular functional pathways as subsets. To identify
the cellular functions regulated by insulin stimulation, we
performed over‐representation analysis of the IRpPs in
the cellular functional pathways. The international protein
index (IPI) (Kersey et al., 2004) IDs of phosphopeptides
in the phosphoproteomic data were converted to KEGG
gene ID using bioDBnet (Mudunuri, Che, Yi, & Stephens,2009) to correspond with the pathway information in the
KEGG database. Over‐representation of the IRpPs for eachpathway was determined using Fisher’s exact test with
FDR using Storey’s procedure (Storey et al., 2004) (significant over‐representation when FDR <0.1). The cellularfunctional pathways in which the IRpPs were significantly
over‐represented were defined as the cellular functional
layer.

4.4 | Inference of protein kinases for protein
phosphorylation
We inferred KSRs for amino acid sequences in the proteins
with quantitatively changed phosphopeptides using a stand‐alone version of NetPhorest (Horn et al., 2014; Miller et al.,2008) (http://netphorest.info/download/netphorest_human.
tsv.xz). NetPhorest inputs are rat protein sequences associated with the IPI in FASTA format (ftp://ftp.ebi.ac.uk/pub/databases/IPI/last_release/current/ipi.RAT.fasta.gz).
The
outputs for NetPhorest are posterior probabilities of an aminoacid residue being recognized by a protein kinase classifier(kinases with similar substrate recognition motifs). Amongthe candidate classifiers, we selected the classifier with thelargest posterior probability value as the kinase classifier

Genes to Cells

  

| 9

related to the amino acid sequence. A predicted KSR is represented as an edge between a kinase classifier as one node
and a phosphorylation site of an IRpP phosphopeptide as theother node. We also extracted individual kinases within ineach class from NetPhorest (view source: http://netphorest.
info/download.shtml) and defined these kinases as responsible protein kinases.

4.5 | Clustering the cellular functional
pathways by the occurrence rates of
responsible protein kinases
We performed hierarchical cluster analysis on the cellularfunctional pathways in the cellular functional layer using theoccurrence rate of the predicted responsible protein kinases.
The occurrence rate of a specific kinase classifier (i) in a specific pathway (j) was calculated as the rate of the numberof quantitatively changed phosphopeptides having the kinaseclassifier i (including the responsible protein kinases) to thetotal number of quantitatively changed phosphopeptides in
the pathway j. The sum of occurrence rates of the kinaseclassifiers in pathway j is 1. We performed the hierarchical clustering of the cellular functions pathways using the
Euclidean distance to calculate the intracluster distances andusing Ward’s method to calculate the intercluster distances(Ward, 1963).

4.6|

Generating motif logos

Motif logos of the quantitatively changed phosphopeptides
included in each cellular functional pathway were generated using enoLOGOS (Workman et al., 2005; http://biodev.hgen.pitt.edu/enologos/) with relative entropy as thelogo plot method. These logos are based on the frequency ofeach amino acid residue within the sequences of the phosphopeptides included in each cellular functional pathway.
Ranges of motifs were provided from −5 to +5 residues
from phosphorylated sites. Statistical tests of amino acidcomposition at each position of the quantitatively changed
phosphopeptides included in the cellular functional layer orin each cluster were performed using iceLogo (Colaert etal., 2009; http://iomics.ugent.be/icelogoserver/index.html)
with percentage difference as the scoring system and a p
value cutoff of 0.05. These logos are based on the statistical significance of frequency of each amino acid residue
within the sequences of the phosphopeptides included ineach cluster against the reference compositions. The Rattus
norvegicus amino acid compositions from Swiss‐Prot andquantitatively changed phosphopeptides included in cellular functions were used as reference compositions for testsof quantitatively changed phosphopeptides included in thecellular functional layer and for tests of those included ineach cluster.

10

|   

Genes to Cells

ACKNOWLEDGMENTS
We thank our laboratory members for critically reading thismanuscript and for their technical assistance with the experiments. The computational analysis of this work was performedin part with support of the supercomputer system of National
Institute of Genetics (NIG) and Research Organization ofInformation and Systems (ROIS). This manuscript was editedby Nancy R. Gough (BioSerendipity, LLC). This work was
supported by the Creation of Fundamental Technologies for
Understanding and Control of Biosystem Dynamics, CREST(JPMJCR12W3), from the Japan Science and Technology
Agency (JST); by the Japan Society for the Promotion of
Science (JSPS) KAKENHI Grant Numbers (17H06300,
17H06299 and 18H03979). K.Y. receives funding from JSPS
KAKENHI Grant Numbers JP15H05582 and JP18H05431,
and “Creation of Innovative Technology for Medical
Applications Based on the Global Analyses and Regulationof Disease‐Related Metabolites,” PRESTO (JPMJPR1538),
from JST. M.F. was supported by JSPS KAKENHI Grant
Number 16K12508. H.K. was supported by JSPS KAKENHI
Grant Number 16H06577. S.U. receives funding from a
Grant‐in‐Aid for Scientific Research on Innovative Areas(18H04801). K.I.N. was supported by the JSPS KAKENHI
Grant Numbers (17H06301 and 17H06306).
CONFLICT OF INTEREST
The authors declare that they have no conflict of interest withthe contents of this article.
DATA AVAILABILIT Y
The dataset used in this analysis including time series ofphosphoproteomic data obtained from Fao cells in response
to insulin stimulation can be obtained from jPOST (accession
number: PXD005900).
ORCID
https://orcid.
Keiichi I. Nakayamaorg/0000-0002-7185-1529
Shinya Kuroda

https://orcid.org/0000-0001-5059-8299

R E F E R E NC E S
Al‐Ayoubi, A. M., Zheng, H., Liu, Y., Bai, T., & Eblen, S. T. (2012).
Mitogen‐activated protein kinase phosphorylation of splicing factor 45 (SPF45) regulates SPF45 alternative splicing site utilization,
proliferation, and cell adhesion. Molecular and Cellular Biology,
32(14), 2880–2893. https://doi.org/10.1128/MCB.06327-11
Alli Shaik, A., Qiu, B., Wee, S., Choi, H., Gunaratne, J., &Tergaonkar, V. (2016). Phosphoprotein network analysis of white

KAWATA et al.

adipose tissues unveils deregulated pathways in response to high‐fat diet. Scientific Reports, 6(1), 25844. https://doi.org/10.1038/
srep25844
Ashburner, M., Ball, C. A., Blake, J. A., Botstein, D., Butler, H.,
Cherry, J. M., … Sherlock, G. (2000). Gene ontology: Tool for theunification of biology. Nature Genetics, 25(1), 25–29. https://doi.
org/10.1038/75556
Boulton, T. G., Nye, S. H., Robbins, D. J., Ip, N. Y., Radziejewska, E.,
Morgenbesser, S. D., … Yancopoulos, G. D. (1991). ERKs: A family of protein‐serine/threonine kinases that are activated and tyrosine phosphorylated in response to insulin and NGF. Cell, 65(4),663–675. https://doi.org/10.1016/0092-8674(91)90098-J
Buescher, J. M., Liebermeister, W., Jules, M., Uhr, M., Muntel, J.,
Botella, E., … Sauer, U. (2012). Global network reorganization during dynamic adaptations of Bacillus subtilis metabolism. Science, 335(6072), 1099–1103. https://doi.org/10.1126/
science.1206871
Chiappino‐Pepe, A., Pandey, V., Ataman, M., & Hatzimanikatis, V.
(2017). Integration of metabolic, regulatory and signaling networks
towards analysis of perturbation and dynamic responses. Current
Opinion in Systems Biology, 2, 59–66. https://doi.org/10.1016/j.
coisb.2017.01.007
Colaert, N., Helsens, K., Martens, L., Vandekerckhove, J., & Gevaert, K.
(2009). Improved visualization of protein consensus sequences by
iceLogo. Nature Methods, 6(11), 786–787. https://doi.org/10.1038/
nmeth1109-786
Croft, D., Mundo, A. F., Haw, R., Milacic, M., Weiser, J., Wu, G.,… D’Eustachio, P. (2014). The reactome pathway knowledgebase. Nucleic Acids Research, 42(D1), D472–D477. https://doi.
org/10.1093/nar/gkt1102
Eblen, S. T. (2018). Extracellular‐regulated kinases: Signaling fromRas to ERK substrates to control biological outcomes. Advancesin Cancer Research, 138, 99–142. https://doi.org/10.1016/
bs.acr.2018.02.004
Emdal, K. B., Pedersen, A.‐K., Bekker‐Jensen, D. B., Tsafou,
K. P., Horn, H., Lindner, S., … Olsen, J. V. (2015). Temporal
proteomics of NGF‐TrkA signaling identifies an inhibitoryrole for the E3 ligase Cbl‐b in neuroblastoma cell differentiation. Science Signaling, 8(374), ra40. https://doi.org/10.1126/
scisignal.2005769
Friedman, A. A., Tucker, G., Singh, R., Yan, D., Vinayagam, A., Hu,
Y., … Perrimon, N. (2011). Proteomic and functional genomic landscape of receptor tyrosine kinase and Ras to extracellular signal‐regulated kinase signaling. Science Signaling, 4(196), rs10. https://doi.
org/10.1126/scisignal.2002029
Giansanti, P., Aye, T. T., van den Toorn, H., Peng, M., van Breukelen,B., & Heck, A. J. R. (2015). An augmented multiple‐protease‐based
human phosphopeptide Atlas. Cell Reports, 11(11), 1834–1843.
https://doi.org/10.1016/j.celrep.2015.05.029
Hartmann, B., Castelo, R., Blanchette, M., Boue, S., Rio, D. C., &Valcárcel, J. (2009). Global analysis of alternative splicing regulation by insulin and wingless signaling in Drosophila cells. Genome
Biology, 10(1), R11. https://doi.org/10.1186/gb-2009-10-1-r11
Hatzimanikatis, V., & Saez‐Rodriguez, J. (2015). Integrative approaches for signalling and metabolic networks. Integrative Biology:
Quantitative Biosciences from Nano to Macro, 7(8), 844–845.
https://doi.org/10.1039/c5ib90030a
Horn, H., Schoof, E. M., Kim, J., Robin, X., Miller, M. L., Diella, F.,… Linding, R. (2014). KinomeXplorer: An integrated platform forKAWATA et al.

kinome biology studies. Nature Methods, 11(6), 603–604. https://doi.org/10.1038/nmeth.2968
Huang, H.‐D., Lee, T.‐Y., Tzeng, S.‐W., & Horng, J.‐T. (2005).
KinasePhos: A web tool for identifying protein kinase‐specificphosphorylation sites. Nucleic Acids Research, 33(Web Serverissue), W226–W229. https://doi.org/10.1093/nar/gki471
Humphrey, S. J., Azimifar, S. B., & Mann, M. (2015). High‐throughput
phosphoproteomics reveals in vivo insulin signaling dynamics. Nature
Biotechnology, 33(9), 990–995. https://doi.org/10.1038/nbt.3327
Humphrey, S. J., James, D. E., & Mann, M. (2015). Protein phosphorylation: A major switch mechanism for metabolic regulation. Trends
in Endocrinology & Metabolism, 26(12), 676–687. https://doi.
org/10.1016/j.tem.2015.09.013
Humphrey, S. J., Yang, G., Yang, P., Fazakerley, D. J., Stöckli, J., Yang,
J. Y., & James, D. E. (2013). Dynamic adipocyte phosphoproteome
reveals that Akt directly regulates mTORC2. Cell Metabolism,
17(6), 1009–1020. https://doi.org/10.1016/j.cmet.2013.04.010
Hyduke, D. R., Lewis, N. E., & Palsson, B. Ø. (2013). Analysis ofomics data with genome‐scale models of metabolism. Molecular
BioSystems, 9(2), 167–174. https://doi.org/10.1039/c2mb25453k
Joyce, A. R., & Palsson, B. Ø. (2006). The model organism as a system: Integrating “omics” data sets. Nature Reviews Molecular Cell
Biology, 7(3), 198–210. https://doi.org/10.1038/nrm1857
Kanehisa, M., Furumichi, M., Tanabe, M., Sato, Y., & Morishima, K.
(2017). KEGG: New perspectives on genomes, pathways, diseases
and drugs. Nucleic Acids Research, 45(D1), D353–D361. https://doi.org/10.1093/nar/gkw1092
Kanehisa, M., Goto, S., Sato, Y., Furumichi, M., & Tanabe, M. (2012).
KEGG for integration and interpretation of large‐scale molecular
data sets. Nucleic Acids Research, 40(Database, issue), D109–D114.
https://doi.org/10.1093/nar/gkr988
Kanshin, E., Kubiniok, P., Thattikota, Y., D’Amours, D., & Thibault,
P. (2015). Phosphoproteome dynamics of Saccharomyces cerevisiae
under heat shock and cold stress. Molecular Systems Biology, 11(6),813. https://doi.org/10.15252/msb.20156170
Kawata, K., Hatano, A., Yugi, K., Kubota, H., Sano, T., Fujii, M., …
Kuroda, S. (2018). Trans‐omic analysis reveals selective responses
to induced and basal insulin across signaling, transcriptional, andmetabolic networks. iScience, 7, 212–229. https://doi.org/10.1016/j.
isci.2018.07.022
Kersey, P. J., Duarte, J., Williams, A., Karavidopoulou, Y., Birney, E., &Apweiler, R. (2004). The International Protein Index: An integrated
database for proteomics experiments. Proteomics, 4(7), 1985–1988.
https://doi.org/10.1002/pmic.200300721
Krüger, M., Kratchmarova, I., Blagoev, B., Tseng, Y.‐H., Kahn, C. R.,& Mann, M. (2008). Dissection of the insulin signaling pathway
via quantitative phosphoproteomics. Proceedings of the National
Academy of Sciences of the United States of America, 105(7), 2451–
2456. https://doi.org/10.1073/pnas.0711713105
Kubota, H., Noguchi, R., Toyoshima, Y., Ozaki, Y., Uda, S., Watanabe,K., … Kuroda, S. (2012). Temporal coding of insulin action throughmultiplexing of the AKT pathway. Molecular Cell, 46(6), 820–832.
https://doi.org/10.1016/j.molcel.2012.04.018
Lange, A. J., Argaud, D., El‐Maghrabi, M. R., Pan, W., Maitra, S. R.,& Pilkis, S. J. (1994). Isolation of a cDNA for the catalytic subunit of rat liver glucose‐6‐phosphatase: Regulation of gene expression in FAO hepatoma cells by insulin, dexamethasone and cAMP.
Biochemical and Biophysical Research Communications, 201(1),302–309. https://doi.org/10.1006/bbrc.1994.1702

Genes to Cells

  

| 11

Lindsay, J. R., McKillop, A. M., Mooney, M. H., Flatt, P. R., Bell, P.
M., & O’Harte, F. P. M. (2003). Meal‐induced 24‐hour profile ofcirculating glycated insulin in type 2 diabetic subjects measured bya novel radioimmunoassay. Metabolism, 52(5), 631–635. https://doi.
org/10.1053/meta.2003.50150
Lizcano, J. M., & Alessi, D. R. (2002). The insulin signalling pathway.
Current Biology: CB, 12(7), R236–R238. https://doi.org/10.1016/
S0960-9822(02)00777-7
Manning, G., Whyte, D. B., Martinez, R., Hunter, T., & Sudarsanam, S.
(2002). The protein kinase complement of the human genome. Science,
298(5600), 1912–1934. https://doi.org/10.1126/science.1075762
Matter, N., Herrlich, P., & König, H. (2002). Signal‐dependent regulation of splicing via phosphorylation of Sam68. Nature, 420(6916),
691–695. https://doi.org/10.1038/nature01153
Miller, M. L., Jensen, L. J., Diella, F., Jorgensen, C., Tinti, M., Li, L., …
Linding, R. (2008). Linear motif atlas for phosphorylation‐dependent signaling. Science Signaling, 1(35), ra2. https://doi.org/10.1126/
scisignal.1159433
Monetti, M., Nagaraj, N., Sharma, K., & Mann, M. (2011). Large‐scale
phosphosite quantification in tissues by a spike‐in SILAC method.
Nature Methods, 8(8), 655–658. https://doi.org/10.1038/nmeth.1647
Mudunuri, U., Che, A., Yi, M., & Stephens, R. M. (2009). bioDBnet:
The biological database network. Bioinformatics (Oxford, England),25(4), 555–556. https://doi.org/10.1093/bioinformatics/btn654
Olsen, J. V., Blagoev, B., Gnad, F., Macek, B., Kumar, C., Mortensen,
P., & Mann, M. (2006). Global, in vivo, and site‐specific phosphorylation dynamics in signaling networks. Cell, 127(3), 635–648.
https://doi.org/10.1016/j.cell.2006.09.026
Olsen, J. V., Vermeulen, M., Santamaria, A., Kumar, C., Miller, M. L.,
Jensen, L. J., … Mann, M. (2010). Quantitative phosphoproteomics reveals widespread full phosphorylation site occupancy duringmitosis. Science Signaling, 3(104), ra3. https://doi.org/10.1126/
scisignal.2000475
Palsson, B., & Zengler, K. (2010). The challenges of integrating multi‐
omic data sets. Nature Chemical Biology, 6(11), 787–789. https://doi.org/10.1038/nchembio.462
Patella, F., Schug, Z. T., Persi, E., Neilson, L. J., Erami, Z., Avanzato,
D., … Zanivan, S. (2015). Proteomics‐based metabolic modelingreveals that fatty acid oxidation (FAO) controls endothelial cell (EC)permeability. Molecular & Cellular Proteomics, 14(3), 621–634.
https://doi.org/10.1074/mcp.M114.045575
Pearson, G., Robinson, F., Beers Gibson, T., Xu, B., Karandikar, M.,
Berman, K., & Cobb, M. H. (2001). Mitogen‐activated protein(MAP) kinase pathways: Regulation and physiological functions.
Endocrine Reviews, 22(2), 153–183. https://doi.org/10.1210/
edrv.22.2.0428
Polonsky, K. S., Given, B. D., & Van Cauter, E. (1988). Twenty‐four‐
hour profiles and pulsatile patterns of insulin secretion in normal
and obese subjects. Journal of Clinical Investigation, 81(2), 442–
448. https://doi.org/10.1172/JCI113339
Reiss, K., Wang, J. Y., Romano, G., Tu, X., Peruzzi, F., & Baserga,
R. (2001). Mechanisms of regulation of cell adhesion and motility
by insulin receptor substrate‐1 in prostate cancer cells. Oncogene,20(4), 490–500. https://doi.org/10.1038/sj.onc.1204112
Sacco, F., Silvestri, A., Posca, D., Pirrò, S., Gherardini, P. F., Castagnoli,
L., … Cesareni, G. (2016). Deep proteomics of breast cancer cells
reveals that metformin rewires signaling networks away from a pro‐
growth state. Cell Systems, 2(3), 159–171. https://doi.org/10.1016/j.
cels.2016.02.005

12

|   

KAWATA et al.

Genes to Cells

Saltiel, A. R., & Kahn, C. R. (2001). Insulin signalling and the regulation of glucose and lipid metabolism. Nature, 414(6865), 799–806.
https://doi.org/10.1038/414799a
Song, C., Ye, M., Liu, Z., Cheng, H., Jiang, X., Han, G., … Zou, H.
(2012). Systematic analysis of protein phosphorylation networks
from phosphoproteomic data. Molecular & Cellular Proteomics,
11(10), 1070–1083. https://doi.org/10.1074/mcp.M111.012625
Storey, J. D., Taylor, J. E., & Siegmund, D. (2004). Strong control, conservative point estimation and simultaneous conservative consistency of false discovery rates: A unified approach. Journal of theRoyal Statistical Society: Series B (Statistical Methodology), 66(1),
187–205. https://doi.org/10.1111/j.1467-9868.2004.00439.x
Tsakiridis, T., Tong, P., Matthews, B., Tsiani, E., Bilan, P. J., Klip,
A., & Downey, G. P. (1999). Role of the actin cytoskeleton in insulin action. Microscopy Research and Technique, 47(2), 79–92.
https://doi.org/10.1002/(SICI)1097-0029(19991015)47:2<79:
AID-JEMT1>3.0.CO;2-S
Vinayagam, A., Kulkarni, M. M., Sopko, R., Sun, X., Hu, Y., Nand,
A., … Perrimon, N. (2016). An integrative analysis of the InR/
PI3K/Akt network identifies the dynamic response to insulin signaling. Cell Reports, 16(11), 3062–3074. https://doi.org/10.1016/j.
celrep.2016.08.029
Ward, J. H. (1963). Hierarchical grouping to optimize an objectivefunction. Journal of the American Statistical Association, 58(301),
236–244. https://doi.org/10.1080/01621459.1963.10500845
Wolf, A., Rietscher, K., Glaß, M., Hüttelmaier, S., Schutkowski, M.,
Ihling, C., … Hatzfeld, M. (2013). Insulin signaling via Akt2
switches plakophilin 1 function from stabilizing cell adhesion topromoting cell proliferation. Journal of Cell Science, 126(Pt 8),
1832–1844. https://doi.org/10.1242/jcs.118992
Workman, C. T., Yin, Y., Corcoran, D. L., Ideker, T., Stormo, G. D.,& Benos, P. V. (2005). enoLOGOS: A versatile web tool for energy normalized sequence logos. Nucleic Acids Research, 33(WebServer, issue), W389–W392. https://doi.org/10.1093/nar/gki439
Yugi, K., Kubota, H., Toyoshima, Y., Noguchi, R., Kawata, K., Komori,Y., … Kuroda, S. (2014). Reconstruction of insulin signal flow from
phosphoproteome and metabolome data. Cell Reports, 8(4), 1171–
1183. https://doi.org/10.1016/j.celrep.2014.07.021
Zhang, Y., Zhang, Y., & Yu, Y. (2017). Global phosphoproteomic
analysis of Insulin/Akt/mTORC1/S6K signaling in rat hepatocytes. Journal of Proteome Research, https://doi.org/10.1021/acs.
jproteome.7b00140

SUPPORTING INFORMATION
Additional supporting information may be found online in
the Supporting Information section at the end of the article.   
How to cite this article: Kawata K, Yugi K, HatanoA, et al. Reconstruction of global regulatory networkfrom signaling to cellular functions using
phosphoproteomic data. Genes Cells. 2018;00:1–12.
https://doi.org/10.1111/gtc.12655