 Research
 Open Access
 Published:
Analysis of significant protein abundance from multiple reactionmonitoring data
BMC Systems Biologyvolume 12, Article number: 123 (2018)
Abstract
Background
Discovering reliable protein biomarkers is one of the most important issues in biomedical research. The ELISA is a traditional technique for accurate quantitation of wellknown proteins. Recently, the multiple reactionmonitoring (MRM) mass spectrometry has been proposed for quantifying newly discovered protein and has become a popular alternative to ELISA. For the MRM data analysis, linear mixed modeling (LMM) has been used to analyze MRM data. MSstats is one of the most widely used tools for MRM data analysis that is based on the LMMs. However, LMMs often provide various significance results, depending on model specification. Sometimes it would be difficult to specify a correct LMM method for the analysis of MRM data. Here, we propose a new logistic regressionbased method for Significance Analysis of Multiple Reaction Monitoring (LRSAM).
Results
Through simulation studies, we demonstrate that LMM methods may not preserve type I error, thus yielding high false positive errors, depending on how random effects are specified. Our simulation study also shows that the LRSAM approach performs similarly well as LMM approaches, in most cases. However, LRSAM performs better than the LMMs, particularly when the effects sizes of peptides from the same protein are heterogeneous. Our proposed method was applied to MRM data for identification of proteins associated with clinical responses of treatment of 115 hepatocellular carcinoma (HCC) patients with the tyrosine kinase inhibitor sorafenib. Of 124 candidate proteins, LMM approaches provided 6 results varying in significance, while LRSAM, by contrast, yielded 18 significant results that were quite reproducibly consistent.
Conclusion
As exemplified by an application to HCC data set, LRSAM more effectively identified proteins associated with clinical responses of treatment than LMM did.
Background
Discovering protein disease biomarkers is an urgently pressing issue in biomedical research [1]. Historically, the enzymelinked immunosorbent assay (ELISA) is a highly accurate protein quantitation technique [2], representing the “gold standard” for measuring levels of specific proteins [3]. However, recent discoveries of many novel proteins, having no available antibodies, now limits the use of the ELISA method [4]. Moreover, for newly discovered proteins, development of high quality ELISA assays requires considerable time and resources [5]. Together, these shortcomings have created a need for a different technique of targeted protein quantitation [6]. One such technique, quantitative mass spectrometry of proteins has advanced significantly over the last decade, and several methods have been developed for relative quantification of targeted proteins, using this technique [7].
In recent years, multiple reaction monitoring (MRM) mass spectrometry has been developed as an attractive tool for targeted proteins and now represents a promising alternative to ELISA for quantification of proteins. To that end, MRM uses sequencespecific tandem mass spectrometer fragmentations of peptides that can provide highly selective measurements for peptidecontaining proteins [8]. Thus, without enrichment or fractionation approaches, MRM assays allow for the quantitation of protein ions within ranges of low (μg/mL) to high (ng/mL) levels [5, 9]. Additionally, development time for MRM assays is relatively shorter and less expensive than that for ELISA, with no requirements (and thus no costs) for antibody development.
Due to these many advantages, MRM assays are being increasingly used in systems biology and clinical investigations [10,11,12,13]. There are some analysis tools to analyze MRM data. Skyline is representative tool to create and analyze MRM dataset [14]. ProteoSign gives some simple statistical analysis results to find differentially expressed proteins [15]. roHitsviz is a webbased tool for visualizing interaction of protein data [16]. However, the development of statistical methods to determine significantly abundant proteins by the MRM assay has not received enough attention, compared to improvement of the MRM assay technology itself [17]. For MRM data analysis, twosample ttests or paired ttests have been applied to identify proteins that change in abundance between two groups [18,19,20]. To test for multiple groups, oneway analysis of variance (ANOVA) has been employed [21, 22].
Recently, a linear mixed model (LMM) approach was proposed for MRM data analysis, and was implemented in MSstats (v3.7.1) [17], resulting in its widespread adoption [23]. LMM approaches include a single parameter for a group effect, peptide effects, and group and peptide interaction effects, along with subject and run effects. LMMs include both subject and run effects that can be designated as either random or fixed, to identify significantly evident proteins by testing the group effect. Recently, we observed that the LMM approach often provides different pvalues for the group effect, from the same data, depending on which effects are treated as random or fixed. In particularly, the LMM test results vary considerably when there is an interaction effect between group and peptide.
In this report, we propose a new logistic regressionbased method for Significance Analysis of Multiple Reaction Monitoring (LRSAM). Unlike LMMs, our LRSAM approach uses a much smaller number of parameters. Moreover, our LRSAM does not require inclusion of all the effects related to the run. Accordingly, our model does not need to specify run effects as random or fixed.
Since LRSAM uses log2transformed relative intensity values, it does not need to include run effects. We consider two separate cases when the effects of peptides in a protein are homogeneous or heterogeneous. For the significance test for proteins, we consider Wald type tests, likelihood ratio tests (LRTs), and score tests [24].
Through various simulation studies, we compared the performance of LMMs to our approach. We first observed cases in which an LMM approach did not preserve type I error, and we also examined cases in which LRSAM methods performed better. Some tests, based on the LRSAM method, were more powerful when the expression pattern of peptide between groups was heterogeneous.
To further establish the translational relevance of our method, we compared LMMs to our approaches in a clinical study to identify candidate serum biomarkers of sorafenib response and prognosis in patients with hepatocellular carcinoma (HCC). Sorafenib is a multikinase inhibitor mainly used for the treatment of various solid tumors including HCC [25]. However, the response rate of sorafenib is about 8% for HCC patients. Due to its high cost, it would be more economical and more beneficial to patients if Sorafenib is applied only to the patients with high chance of response. Moreover, there are no good markers to predict the patients’ response to sorafenib [26]. From May 2013 to August 2014, 115 HCC patient serum samples were collected from Seoul National University Hospital as part of an ongoing HCC study. One hundred twentyfour candidate protein biomarkers, and hepatic diseaseassociated proteins, were chosen from 50,265 proteins, based on the LiverAtlas database [27]. Sorafenib responses were measured according to the modified response evaluation criteria in solid tumors (mRECIST) guidelines [28]. Patients with complete responses, partial responses, or stable disease were categorized as responders, while those with progressive disease were categorized nonresponders. Both LMM approaches and LRSAM methods were employed to analyze the MRM data.
Methods
MRM dataset structure
Table 1 shows the structure of MRM dataset. Sample ID indicates individual, Run the shared run membership of the endogenous, and Area Ratio the expression of peptide detected. The MRM data has a hierarchical structure such that one protein contains several peptides. Since our goal is to identify proteins that are significantly different between the group, we need to combine the summarized expression of peptides from the same protein efficiently. Also, the batch effect would occur when MRM experiments were performed in different batches separately. These characteristics of MRM make it difficult to use the standard analyses such as two sample ttest.
Review of LMM approach
For a given protein, suppose it is comprised of K peptides. Let y_{i, j(i), k, l} denote the log2(intensity) value of the jth subject, nested in the ith group of the kth peptide and the lth run. Then, the LMM used in MSstats is given as follows:
where μ is the global mean; G_{i} is the ith group effect; S(G)_{j(i)} representing the jth subject effect nested in the ith group; P_{k} stands for the kth peptide effect; R_{l} stands for the lth run effect, (G × P)_{ik} stands for the interaction effect between the ith group and the kthe peptide; and (P × R)_{kl} signifies the interaction effect between the kth peptide and the lth run. When all the effects are treated as fixed, these parameters have the following restrictions: \( \sum \limits_{i=0}^2{G}_i=0 \), \( \sum \limits_{j(i)=1}^{J(i)}S{(G)}_{j(i)}=0 \), \( \sum \limits_{k=1}^K{P}_k=0 \), \( \sum \limits_{l=1}^L{R}_l=0 \), \( \sum \limits_{i=0}^2{\left(G\times P\right)}_{i,k}=0 \), \( \sum \limits_{k=1}^K{\left(G\times P\right)}_{i,k}=0 \), \( \sum \limits_{k=1}^K{\left(P\times R\right)}_{k,l}=0,\sum \limits_{l=1}^L{\left(P\times R\right)}_{k,l}=0 \), and \( {\epsilon}_{i,j(i),k,l}\sim N\left(0,{\sigma}_{\epsilon}^2\right) \). Here, G_{0} stands for the effect of a reference group of MRM data. When the subject and run effects are treated as random, the restrictions on S(G)_{j(i)}, R_{l}, and (P × R)_{kl} are replaced by \( S{(G)}_{j(i)}\sim N\left(0,{\sigma}_S^2\right) \), \( {R}_l\sim N\left(0,{\sigma}_R^2\right) \) and \( {\left(P\times R\right)}_{kl}\sim N\left(0,{\sigma}_{P\times R}^2\right) \), respectively.
For most MRM data analyses, the investigator’s interest lies in determining which proteins differ in abundance between two groups. Thus, the hypothesis of interest is given below for comparing two groups:
The MSstats uses the ttest for this hypothesis.
LRSAM approach
Our proposed LRSAM approach uses a log2transformed relative intensity value instead of the original log2transformed intensity value itself. A log2transformed relative intensity value was derived as y_{i, j(i), k, l} − y_{0, 0(0), k, l}. In that scenario, P_{k}, R_{l} and (P × R)_{kl} effects that share the same k and l values are removed. If we denote y_{j, k} as a log2transformed relative intensity value of the jth subject and the kth peptide, where j = 1,J(1),J(1) + 1,…,J(1) + J(2). Then, the log2transformed relative intensity value y_{j, k} stands for the y_{1, j(1), k, l} − y_{0, 0(0), k, l} for j = 1,…,J(1), and it stands for the y_{2, j(2), k, l} − y_{0, 0(0), k, l} for j = J(1) + 1,…,J(1) + J(2).
LRSAM with fixed effect
Consider a logistic regression model using log2transformed relative intensity value, as follows:
where Z_{j} is a group indicator of the jth subject that is assumed to follow a Bernoulli distribution; α is an intercept; and β_{k} is the coefficient of the kth peptide. Note that the β_{k} values are related to G_{1} − G_{2} + (G × P)_{1, k} − (G × P)_{2, k} of model (1). Therefore, if we treat all β_{k}’s as fixed, the hypothesis for comparing the two groups is:
To test (4), we consider the likelihood ratio test (L) with K degrees of freedom, as follows:
Here, \( {\widehat{\alpha}}_0 \) is the maximum likelihood estimate (MLE) of α under the null hypothesis (see above). \( {l}_0\left({\widehat{\alpha}}_0\right) \) is the maximum likelihood value under the null hypothesis, and \( \widehat{\alpha},{\widehat{\beta}}_1,\dots, {\widehat{\beta}}_K \) are the MLEs of α, β_{1}, …, β_{K}, respectively. It is also known that L asymptotically follows a chisquare distribution, with K degrees of freedom.
A Wald type test (W) statistic for analysis, is given below:
Here, W also asymptotically follows a chisquare distribution, with K degrees of freedom under the null hypothesis. If we assume that the β_{k} values are homogeneous, a Wald test with 1 degree of freedom (W_{1}) can be considered as follows:
in which \( {\widehat{\beta}}_p \) is the weighted average of \( {\widehat{\beta}}_k \) as \( \sum \limits_{k=1}^K{t}_k{\widehat{\beta}}_k \), where t_{k} is a weight given as \( \frac{1/ Var\left({\widehat{\beta}}_k\right)}{\sum_{k=1}^K1/ Var\left({\widehat{\beta}}_k\right)} \), and A is the variance of \( {\widehat{\beta}}_p \), given as \( {\left(\begin{array}{c}{t}_1\\ {}\vdots \\ {}{t}_K\end{array}\right)}^{\prime } Var\left({\widehat{\beta}}_1,\dots, {\widehat{\beta}}_K\right)\left(\begin{array}{c}{t}_1\\ {}\vdots \\ {}{t}_K\end{array}\right) \). Thus, W_{1} asymptotically follows a chisquare distribution, with 1 degree of freedom, under the null hypothesis. Moreover, if the β_{k}values are homogeneous, Eq. (3) (shown above) can be reduced to the following model:
The hypothesis of interest from this model (8) is H_{0} : β^{∗} = 0, and a Wald test (W_{S}) is given by:
Here, \( {\widehat{\beta}}^{\ast } \) is the maximum likelihood estimate of β^{∗} from the model (8), and W_{S} asymptotically follows a chisquare distribution, with 1 degree of freedom.
LRSAM with random effects
When effects of peptides were heterogeneous, we could detect significant heterogenous effect by assuming random effects on coefficients β_{k}. This assumption is commonly used in metaanalysis or rare variant analysis of genetic study [24, 29]. If we assume that β_{k} follows a normal distribution, with mean 0 and variance w_{k}τ (where w_{k} is a known prior weight), then the logistic regression model (3) is expanded to a mixed effect model, and the hypothesis (4) is equivalent to:
Since the hypothesis H_{0} : τ = 0 is on the boundary of the parameter space, the variancecomponent score test can be considered. The score test statistic of the variancecomponent for (10) is:
Here, \( {\widehat{\boldsymbol{\mu}}}_0 \) is the estimated probability under H_{0}; K = YWY^{′}, where Y = [Y_{1}, …, Y_{K}] and Y_{k} = (y_{1k}, …, y_{nk})^{′}; Z means the group indicator vector, and W is a diagonal matrix with the kth element as w_{k}.
It is known that S_{VC} follows a mixture of chisquare distributions \( {\sum}_{k=1}^K{\uplambda}_k{\chi_{1,k}}^2 \), where χ_{1, k}^{2}’s are independent chisquare distributions with 1 degree of freedom, and λ_{k} is the kth eigenvalue of P^{1/2}KP^{1/2} [24]. Here, \( \boldsymbol{P}={\widehat{\boldsymbol{V}}}^{\mathbf{1}}{\widehat{\boldsymbol{V}}}^{\mathbf{1}}\mathbf{1}\left({\mathbf{1}}^{\prime }{\widehat{\boldsymbol{V}}}^{\mathbf{1}}\mathbf{1}\right){\mathbf{1}}^{\prime }{\widehat{\boldsymbol{V}}}^{\mathbf{1}} \), where \( \widehat{\boldsymbol{V}} \) is a diagonal matrix with the kth element as \( {\widehat{\mu}}_{0k}\left(1{\widehat{\mu}}_{0k}\right) \). For simplicity, we assume a flat prior weight given by w_{k} = 1 for k = 1, …, K.
Simulation design
We next performed simulation studies to investigate the power of the LMM and LRSAM approaches for comparing two groups, and whether they preserve type I error. For analysis, there were four LMMs available, depending on how the random or fixed effects were specified: (i) LMM(FF), with fixed subject and run effects, (ii) LMM(FR), with fixed subject effect and random run effects, (iii) LMM(RF), with random subject and fixed run effects, and (iv) LMM(RR), with random subject and run effects. For each simulated dataset, the best LMM, LMM(best), was selected among four LMMs with the smallest Akaike information criterion (AIC) value [30]. Thus, there were five LRSAM test statistics, L, W, W_{1}, W_{S}, and S_{VC}, for comparison.
Simulation data was generated from the model, using 1000 repetitions. The global mean, μ, was arbitrarily set to 15, while the reference effects, G_{0}, S(G)_{0, (0)}, and (G × P)_{0, k}, for k = 1, …, K, were set to 0. The normal distribution, with mean 0 and variance 0.5, was set as the error distribution, and the number of peptides was set to 4. For the random subject effect, we generated S(G)_{j(i)} from the identical normal distribution independently, with mean 0 and variance 0.25 for i = 1, 2. For the fixed subject effect, we set S(G)_{j(i)} as follows:
The peptide effect, P_{k}, was considered as a fixed effect, and was set as follows:
For the random run effect, we generated R_{l} from the identical normal distribution, independently, with mean 0 and variance 0.25. For the fixed run effect, we set R_{l} as follows:
\( {R}_l={e}_R+\frac{2\left(m1\right)}{M1}{e}_R\kern0.5em \mathrm{for}\kern0.5em m=1,\dots, M,\kern0.5em \mathrm{where}\kern0.5em {e}_R>0\kern0.5em \mathrm{and}\kern0.5em {e}_R=\sqrt{3\frac{M1}{M+1}}{\sigma}_R,\kern0.5em \mathrm{with}\kern0.5em {\sigma}_R^2=0.25 \)
We then considered the peptide×run interaction effect as either random or fixed, followed by the type of run effect. For the random peptide by run interaction effect, we generated (P × R)_{k, l} from the identical normal distribution independently, with mean 0 and variance 0.1. For the fixed peptide×run interaction effect, we set (P × R)_{k, l} as follows:
Provided those parameters, we further considered four group effect scenarios (GSs), as follows:
For detecting the group×peptide interaction effect, we set (G × P)_{i, k}, as follows:
Here, using σ_{GP} values determined from the squared average of the interaction effect, we considered three interaction effect scenarios (ISs) for σ_{GP}, as follows:
Materials
To identify candidate serum hepatocellular carcinoma (HCC) biomarkers for prognosis and response to the tyrosine kinase inhibitor sorafenib, data from 115 HCC patients who had undergone continuous administration of sorafenib for more than 6 weeks, were collected from Seoul National University Hospital, as part of an ongoing study between May 2013 and August 2014. HCC was diagnosed by histological or radiological evaluation, with reference to the American Association for the Study of Liver Diseases (AASLD) [31] or the European Association for the Study of the Liver (EASL) [32] guidelines. All procedures/analyses were approved by the Seoul National University Hospital Institutional Review Board (IRB protocol No. 0506–150005). Sorafenib response was evaluated using the modified response evaluation criteria in solid tumors (mRECIST) [28], using independent radiologic assessments. Patients with complete response, partial response, and stable disease were categorized as responders, while those with progressive disease were categorized as nonresponders.
Toward these objectives, a total of 115 serum samples were randomly separated into two batches: the 1st batch consisted of 65 samples and the 2nd batch consisted of 50 samples. Of those, the total number of responders was 40, and the number of nonresponders was 75. One hundred twentyfour candidate protein biomarkers, known hepatic diseaseassociated proteins, were chosen from 50,265 proteins, based on the LiverAtlas database [27]. One to seven peptides, comprised within each protein, were measured, and counted (Table 2).
Since there were some cases in which LMM(FF), LMM(FR), and L did not preserve type I error in our simulation studies, we applied LMM(RF), LMM(RR), W, W_{1}, W_{S} and S_{VC} to analyze the MRM data. Quantile normalization, provided within the MSstats package [20], was employed for preprocessing. To adjust the mean differences between the two batches, batch indicators were included in models (1), (3), and (8) for MRM data analysis.
Results
Simulation results
The type I error and the empirical power were estimated as the proportion of pvalues under 0.05, out of 1000 repetitions. The type I error rates of LMMs and LRSAM are summarized in Table 3. The true type of subject effect more strongly affected the type I error rate of LMMs, as compared to the run effect. Analogously, when the true subject effect was fixed, the four LMMs controlled the type I errors. Among the five LRSAM tests, L could not control the type I error rate when the sample size was 20, while the other four LRSAM tests could. When the sample sizes were 50 and 100, all LRSAM tests controlled type I error. When the true subject effect was random, LMM(FF) and LMM(FR) could not control type I error, whereas LMM(RF) and LMM(RR) could. L did not control the type I error rate when the sample sizes were 20 and 50, while the other four LRSAM tests did. The AIC value of LMM(FF) tended to be the smallest among the four LMMs, under any conditions. Thus, LMM(FF) was most frequently selected as the best LMM.
Some LRSAM tests showed increased power as the interaction effects became large. When there was only an interaction effect, without group effects, the power of the LMMs W_{1}, and W_{S} did not increase, as shown in Table 4. As the sample size and the interaction effect became large, W and S_{VC} showed increased power. Moreover, S_{VC} provided higher power than W, when the sample size was 50, while W provided higher power than S_{VC}, when the sample size was 100.
The size of the group effect showed large effects on all LMMs and LRSAM methods. The estimated powers of the five LMMs, and the LRSAM methods for GS 1, are shown in Table 5. The powers of L, W, W_{1}, and S_{VC} were all affected by the size of the interaction effect, when the group effect scenario was GS 1, while those of the five LMMs did not. As the interaction effect increased, the powers of L, W, and S_{VC} increased, while the powers of four of the LMMs, and W_{S}, did not, and the power of the W_{1} method even decreased. When the interaction scenario was IS 1, the powers of the L and S_{VC} methods were higher than those of the LMMs, for a sample size of 100. For the IS 2 scenario, once the sample size exceeded 50, the L and S_{VC} methods provided higher powers than those of the LMMs.
Collectively, as the size of group effects became large, all methods provided increased power. As the interaction effect became large, only the L, W, and S_{VC} methods provided increased power.
HCC data analysis
Using the approaches described above, the maximum –log10transformed pvalues of the linear mixed models (LMMs) were higher than those of the LRSAM methods. However, the correlation of the transformed pvalues of LMM(RF) and LMM(RR) was lower than that of LRSAM methods, as shown in Fig. 1. The correlation between LMM(RF) and LMM(RR) was 0.6558, while the minimum correlation among the LRSAM methods was 0.784. The correlation between W_{S} and S_{VC} was the highest (0.9627), among the LRSAM methods.
The lowest pvalues of LMM(RF) and LMM(RR) were lower than those of the LRSAM methods, valued at 1.87 × 10^{−18} and 1.09 × 10^{−30}, respectively. The lowest pvalue of W, W_{1}, and W_{S} was 2 × 10^{−7}, and the lowest pvalue of S_{VC} was 3.49 × 10^{−14} (Table 6). Although the lowest pvalues of the LRSAM methods were less significant than those of the LMMs, and for many proteins, LRSAM methods provided higher significance results (Fig. 1).
The LMM(RF) and LMM(RR) methods identified 11 and 37 proteins, respectively. There were six (14.29%) proteins simultaneously identified by LMM(RF) and LMM(RR). For LRSAM methods, W, W_{1}, W_{S}, and S_{VC} identified 28, 19, 22, and 29 proteins, respectively. Among these proteins 18 (52.94%) were simultaneously identified by four methods. Finally, there were four (7.84%) proteins identified by all six methods of LMM and LRSAM (Fig. 2).
All LMMs and LRSAM methods simultaneously provided significance results for the presence of the proteins GPX3, IGHG1, IGHG3 and IGJ. Among these, IGJ (immunoglobulin J chain, linker protein for immunoglobulin alpha and mu) was previously reported as having a significant difference of expression in HCC tumors [33], while GPX3 (glutathione peroxidase 3) was reported as a tumor suppressor in HCC [34].
There were 14 proteins that the LRSAM methods, but not the LMMs, simultaneously provided significance results of expression (Table 6). Among those 14, seven were previously reported, including FBLN1 (Fibulin1), a tumor suppressor gene in HCC [35], SHBG (sex hormonebinding globulin), a prediagnostic risk marker for HCC [36], and LG3BP (galectin3binding protein), a potential marker in six cancer types, including HCC, lymphoma, NPC (nasopharyngeal carcinoma), CRC (colorectal carcinoma), and oral cancers [37]. Additionally, CATB (Cathepsin B) was previously reported as a potential candidate cancer biomarker in HCC [38], HPT (haptoglobin) was reported to associate with tumor progression in HCC [39], while POSTN (periostin) was reported as a marker for malignant transformation of hepatocytes [40]. Similarly, 14–33S (14–33 protein sigma) was reported as downregulated in HCC [41].
There were two proteins that LMMs simultaneously provided a significance result, while LRSAM methods did not. Resultant pvalues for the simultaneously identified 20 proteins are shown in Table 6. Two proteins that demonstrated significance by LMM(RF) and LMM(RR) were also previously reported: CD5L (CD5 antigenlike), reported as differentially expressed in hepatitis C patients [42]: and APOA4, (apolipoprotein AIV), reported as misexpressed in liver metabolic disorders [43].
Discussion
We examined the performance of LMMs and LRSAM methods, through extensive simulation studies. When the true subject effect was random, LMM(FF) and LMM(FR) did not preserve type I error (Table 3). We also applied Akaike information criterion (AIC), for model selection, to check whether or not the bestperforming LMM preserved type I error. However, our empirical study showed that LMM(FF) had the smallest AIC value for most simulation settings, which made it difficult to use AIC as a model selection criterion. On the other hand, LRSAM methods, except L, well preserved type I error under any type of subject and run effects.
Since the hypothesis (2) did not consider the interaction effect, the power of the LMM approach was not affected under any size of the interaction effect. However, when the interaction effect (without the group effect) was considered, hypothesis (2) could not perform well, as we observed through our simulation studies. On the other hand, testing hypotheses (8) with W and (10) with S_{VC}, changes under any size group effect and interaction effect were detectable. Additionally, LRSAM methods did not provide higher power than LMMs, when there were no or only weak interaction effects.
Although LMM provided more protein significance results than LRSAM methods, the proteins identified by the LMM(RF) and LMM(RR) approaches were very different from each other. On the other hand, LRSAM methods provided more consistent significance results than those of the LMMs. Moreover, LRSAM methods provided greater significance results than those identified by LMM(RF), for most proteins, and several proteins were identified only by LRSAM methods. However, there was a still performance difference between the choices of fixed and random effect models in LRSAM methods, as was in the LMM. This difference was caused by different testing hypotheses between W and S_{VC} in LRSAM. In addition, note that S_{VC} allows the heterogenous effects of peptides, while W does not.
Conclusion
LMMs have been widely used to identify proteins with significantly altered abundance, in distinct disease states, based on MRM assays [39]. However, we found that LMM approaches provide inconsistent significance results for the same MRM data, depending on which effects are treated as random or fixed by simulation results. It is a wellknown property of LMMs that the variance of model parameters are underestimated, when the fixed effect model is fitted but the true effects are random, and vice versa [44]. As a result, the protein significance results of LMMs may vary, depending on whether the true effect is random or fixed, and we also observed this phenomenon, as shown in Fig. 1 and Table 6. Thus, it is highly important to correctly specify the effect as random or fixed.
In this paper, we propose a new logistic regressionbased method for Significance Analysis of Multiple Reaction Monitoring (LRSAM). Unlike LMMs, our LRSAM approach uses a much smaller number of parameters. Moreover, our LRSAM does not require inclusion of all the effects related to the run. Accordingly, our model does not need to specify run effects as random or fixed.
In simulation study, our LRSAM preserved type I error when the true subject effect was random, while LMM(FF) and LMM(FR) did not. In real data analysis, although LMM provided more protein significance results than LRSAM methods, LRSAM methods provided more consistent significance results than those of the LMMs. Thus, our proposed method, LRSAM could give more reliable results than previous protein studies.
Abbreviations
 AASLD:

American Association for the Study of Liver Diseases
 AIC:

Akaike information criterion
 ANOVA:

Analysis of variance
 EASL:

European Association for the Study of the Liver
 ELISA:

Enzymelinked immunosorbent assay
 GS:

Group effect scenarios
 HCC:

Hepatocellular carcinoma
 IS:

Interaction effect scenarios
 LMM:

Linear mixed model
 mRECIST:

Modified response evaluation criteria in solid tumors
 MRM:

Multiple reaction monitoring
References
 1.
Frantzi M, Bhat A, Latosinska A. Clinical proteomic biomarkers: relevant issues on study design & technical considerations in biomarker development. Clin Transl Med. 2014;3(1):7.
 2.
Shi T, Su D, Liu T, Tang K, Camp DG, Qian WJ, Smith RD. Advancing the sensitivity of selected reaction monitoringbased targeted quantitative proteomics. Proteomics. 2012;12(8):1074–92.
 3.
Pan S, Aebersold R, Chen R, Rush J, Goodlett DR, McIntosh MW, Zhang J, Brentnall TA. Mass spectrometry based targeted protein quantification: methods and applications. J Proteome Res. 2008;8(2):787–97.
 4.
Shi T, Fillmore TL, Sun X, Zhao R, Schepmoes AA, Hossain M, Xie F, Wu S, Kim JS, Jones N. Antibodyfree, targeted massspectrometric approach for quantification of proteins at low picogram per milliliter levels in human plasma/serum. Proc Natl Acad Sci. 2012;109(38):15395–400.
 5.
Haab BB, Paulovich AG, Anderson NL, Clark AM, Downing GJ, Hermjakob H, LaBaer J, Uhlen M. A reagent resource to identify proteins and peptides of interest for the Cancer community a WORKSHOP REPORT. Mol Cell Proteomics. 2006;5(10):1996–2007.
 6.
Rifai N, Gillette MA, Carr SA. Protein biomarker discovery and validation: the long and uncertain path to clinical utility. Nat Biotechnol. 2006;24(8):971–83.
 7.
Mesri M. Advances in proteomic technologies and its contribution to the field of Cancer. Adv Med. 2014;2014:238045.
 8.
Liebler DC, Zimmerman LJ. Targeted quantitation of proteins by mass spectrometry. Biochemistry. 2013;52(22):3797–806.
 9.
Kuzyk MA, Smith D, Yang J, Cross TJ, Jackson AM, Hardie DB, Anderson NL, Borchers CH. Multiple reaction monitoringbased, multiplexed, absolute quantitation of 45 proteins in human plasma. Mol Cell Proteomics. 2009;8(8):1860–77.
 10.
Hale JE. Advantageous uses of mass spectrometry for the quantification of proteins. Int J Proteomics. 2013;2013:219452.
 11.
Carr SA, Abbatiello SE, Ackermann BL, Borchers C, Domon B, Deutsch EW, Grant RP, Hoofnagle AN, Hüttenhain R, Koomen JM. Targeted peptide measurements in biology and medicine: best practices for mass spectrometrybased assay development using a fitforpurpose approach. Mol Cell Proteomics. 2014;13(3):907–17.
 12.
Picotti P, Aebersold R. Selected reaction monitoringbased proteomics: workflows, potential, pitfalls and future directions. Nat Methods. 2012;9(6):555–66.
 13.
Grebe SK, Singh RJ. LCMS/MS in the clinical laboratory–where to from here? Clin Biochem Rev. 2011;32(1):5.
 14.
MacLean B, Tomazela DM, Shulman N, Chambers M, Finney GL, Frewen B, Kern R, Tabb DL, Liebler DC, MacCoss MJ. Skyline: an open source document editor for creating and analyzing targeted proteomics experiments. Bioinformatics. 2010;26(7):966–8.
 15.
Efstathiou G, Antonakis AN, Pavlopoulos GA, Theodosiou T, Divanach P, Trudgian DC, Thomas B, Papanikolaou N, Aivaliotis M, Acuto O, et al. ProteoSign: an enduser online differential proteomics statistical analysis platform. Nucleic Acids Res. 2017;45(W1):W300–6.
 16.
Knight JDR, Choi H, Gupta GD, Pelletier L, Raught B, Nesvizhskii AI, Gingras AC. ProHitsviz: a suite of web tools for visualizing interaction proteomics data. Nat Methods. 2017;14(7):645–6.
 17.
Chang CY, Picotti P, Hüttenhain R, HeinzelmannSchwarz V, Jovanovic M, Aebersold R, Vitek O. Protein significance analysis in selected reaction monitoring (SRM) measurements. Mol Cell Proteomics. 2012;11(4):M111. 014662.
 18.
Huang SM, Bisogno T, Trevisani M, AlHayani A, De Petrocellis L, Fezza F, Tognetto M, Petros TJ, Krey JF, Chu CJ. An endogenous capsaicinlike substance with high potency at recombinant and native vanilloid VR1 receptors. Proc Natl Acad Sci. 2002;99(12):8400–5.
 19.
Zemann B, Kinzel B, Müller M, Reuschel R, Mechtcheriakova D, Urtz N, Bornancin F, Baumruker T, Billich A. Sphingosine kinase type 2 is essential for lymphopenia induced by the immunomodulatory drug FTY720. Blood. 2006;107(4):1454–8.
 20.
Mani D, Abbatiello SE, Carr SA. Statistical characterization of multiplereaction monitoring mass spectrometry (MRMMS) assays for quantitative proteomics. BMC Bioinformatics. 2012;13(Suppl 16):S9.
 21.
Yassine HN, Jackson AM, Reaven PD, Nedelkov D, Nelson RW, Lau SS, Borchers CH. The application of multiple reaction monitoring to assess ApoAI methionine oxidations in diabetes and cardiovascular disease. Transl Proteomics. 2014;4:18–24.
 22.
Zhang P, Kirk JA, Ji W, dos Remedios CG, Kass DA, Van Eyk JE, Murphy AM. Multiple reaction monitoring to identify sitespecific troponin I phosphorylated residues in the failing human heart. Circulation. 2012. https://doi.org/10.1161/CIRCULATIONAHA.112.096388.
 23.
Pursiheimo A, Vehmas AP, Afzal S, Suomi T, Chand T, Strauss L, Poutanen M, Rokka A, Corthals GL, Elo LL. Optimization of statistical methods impact on quantitative proteomics data. J Proteome Res. 2015;14(10):4118–26.
 24.
Wu MC, Lee S, Cai T, Li Y, Boehnke M, Lin X. Rarevariant association testing for sequencing data with the sequence kernel association test. Am J Hum Genet. 2011;89(1):82–93.
 25.
Ben Mousa A. Sorafenib in the treatment of advanced hepatocellular carcinoma. Saudi J Gastroenterol. 2008;14(1):40–2.
 26.
Kim H, Yu SJ, Yeo I, Cho YY, Lee DH, Cho Y, Cho EJ, Lee JH, Kim YJ, Lee S, et al. Prediction of response to Sorafenib in hepatocellular carcinoma: a putative marker panel by multiple reaction monitoringmass spectrometry (MRMMS). Mol Cell Proteomics. 2017;16(7):1312–23.
 27.
Zhang Y, Yang C, Wang S, Chen T, Li M, Wang X, Li D, Wang K, Ma J, Wu S. LiverAtlas: a unique integrated knowledge database for systemslevel research of liver and hepatic disease. Liver Int. 2013;33(8):1239–48.
 28.
Camacho JC, Kokabi N, Xing M, Prajapati HJ, ElRayes B, Kim HS. Modified response evaluation criteria in solid tumors and European Association for the Study of the liver criteria using delayedphase imaging at an early time point predict survival in patients with unresectable intrahepatic cholangiocarcinoma following yttrium90 radioembolization. J Vasc Interv Radiol. 2014;25(2):256–65.
 29.
Riley RD, Higgins JP, Deeks JJ. Interpretation of random effects metaanalyses. BMJ. 2011;342:d549.
 30.
Akaike H. Information theory and an extension of the maximum likelihood principle. In: Selected papers of Hirotugu Akaike. New York: Springer; 1998. p. 199–213.
 31.
Bruix J, Sherman M. Management of hepatocellular carcinoma: an update. Hepatology. 2011;53(3):1020–2.
 32.
Ronot M, Bouattour M, Wassermann J, Bruno O, Dreyer C, Larroque B, Castera L, Vilgrain V, Belghiti J, Raymond E. Alternative response criteria (Choi, European association for the study of the liver, and modified response evaluation criteria in solid tumors [RECIST]) versus RECIST 1.1 in patients with advanced hepatocellular carcinoma treated with sorafenib. Oncologist. 2014;19(4):394–402.
 33.
Neo SY, Leow CK, Vega VB, Long PM, Islam AF, Lai P, Liu ET, Ren EC. Identification of discriminators of hepatoma by gene expression profiling using a minimal dataset approach. Hepatology. 2004;39(4):944–53.
 34.
Qi X, Ng K, Lian QZ, Liu XB, Li CX, Geng W, Ling CC, Ma YY, Yeung WH, Tu WW. Clinical significance and therapeutic value of glutathione peroxidase 3 (GPx3) in hepatocellular carcinoma. Oncotarget. 2014;5(22):11103–20.
 35.
Kanda M, Nomoto S, Okamura Y, Hayashi M, Hishida M, Fujii T, Nishikawa Y, Sugimoto H, Takeda S, Nakao A. Promoter hypermethylation of fibulin 1 gene is associated with tumor progression in hepatocellular carcinoma. Mol Carcinog. 2011;50(8):571–9.
 36.
Lukanova A, Becker S, Hüsing A, Schock H, Fedirko V, Trepo E, Trichopoulou A, Bamia C, Lagiou P, Benetou V. Prediagnostic plasma testosterone, sex hormonebinding globulin, IGFI and hepatocellular carcinoma: etiological factors or risk markers? Int J Cancer. 2014;134(1):164–73.
 37.
Wu CC, Hsu CW, Chen CD, Yu CJ, Chang KP, Tai DI, Liu HP, Su WH, Chang YS, Yu JS. Candidate serological biomarkers for cancer identified from the secretomes of 23 cancer cell lines and the human protein atlas. Mol Cell Proteomics. 2010;9(6):1100–17.
 38.
Luk JM, Lam BY, Lee NP, Ho DW, Sham PC, Chen L, Peng J, Leng X, Day PJ, Fan ST. Artificial neural networks and decision tree model analysis of liver cancer proteomes. Biochem Biophys Res Commun. 2007;361(1):68–73.
 39.
Ang IL, Poon TC, Lai PB, Chan AT, Ngai SM, Hui AY, Johnson PJ, Sung JJ. Study of serum haptoglobin and its glycoforms in the diagnosis of hepatocellular carcinoma: a glycoproteomic approach. J Proteome Res. 2006;5(10):2691–700.
 40.
Riener MO, Fritzsche FR, Soll C, Pestalozzi BC, ProbstHensch N, Clavien PA, Jochum W, Soltermann A, Moch H, Kristiansen G. Expression of the extracellular matrix protein periostin in liver tumours and bile duct carcinomas. Histopathology. 2010;56(5):600–6.
 41.
Wang X, Chen Y, Qb H, Cy C, Wang H, Liu Z, CHk C, Yew DT, Lin M, He M. Proteomic identification of molecular targets of gambogic acid: role of stathmin in hepatocellular carcinoma. Proteomics. 2009;9(2):242–53.
 42.
Gangadharan B, Antrobus R, Dwek RA, Zitzmann N. Novel serum biomarker candidates for liver fibrosis in hepatitis C patients. Clin Chem. 2007;53(10):1792–9.
 43.
Graveel CR, Jatkoe T, Madore SJ, Holt AL, Farnham PJ. Expression profiling and identification of novel genes in hepatocellular carcinomas. Oncogene. 2001;20(21):2704–12.
 44.
Chung Y, RabeHesketh S, Choi IH. Avoiding zero betweenstudy variance estimates in randomeffects metaanalysis. Stat Med. 2013;32(23):4071–89.
Acknowledgements
Not applicable
Funding
This research was supported by a grant of the Korea Health Technology R&D Project through the Korea Health Industry Development Institute (KHIDI), funded by the Ministry of Health & Welfare, Republic of Korea (grant number: HI16C2037, HI15C2165). Publication of this article was sponsored by HI16C2037 grant.
Availability of data and materials
Not applicable.
About this supplement
This article has been published as part of BMC Systems Biology Volume 12 Supplement 9, 2018: Proceedings of the 29th International Conference on Genome Informatics (GIW 2018): systems biology. The full contents of the supplement are available online at https://bmcsystbiol.biomedcentral.com/articles/supplements/volume12supplement9.
Author information
Affiliations
Contributions
JJ and TP developed statistical method. JJ performed statistical analysis. HK, SJY, EJC, JH Y, YSK and TP conceived and planned the experiments. HK, JP, and JJY carried out the experiments. JJ, JG. YK and TP planned and carried out the simulations. YYC, DHL, YC, EJC, JHL, and YJK contributed to sample preparation. JJ, JG, YK, YSK, TP contributed to the interpretation of the results. JJ and TP took the lead in writing the manuscript. All authors have read and approved the final manuscript.
Corresponding author
Correspondence to Taesung Park.
Ethics declarations
Ethics approval and consent to participate
To identify candidate serum hepatocellular carcinoma (HCC) biomarkers for prognosis and response to the tyrosine kinase inhibitor sorafenib, data from 115 HCC patients who had undergone continuous administration of sorafenib for more than 6 weeks, were collected from Seoul National University Hospital, as part of an ongoing study between May 2013 and August 2014. HCC was diagnosed by histological or radiological evaluation, with reference to the American Association for the Study of Liver Diseases (AASLD) [31] or the European Association for the Study of the Liver (EASL) [32] guidelines. All procedures/analyses were approved by the Seoul National University Hospital Institutional Review Board (IRB protocol No. 0506–150005).
Consent for publication
Not applicable.
Competing interests
The authors declare that they have no competing interests.
Publisher’s Note
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Rights and permissions
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
About this article
Published
DOI
Keywords
 Multiple reactionmonitoring (MRM)
 Protein
 Logistic regressionbased method for significance analysis of multiple reaction monitoring (LRSAM)
 Hepatocellular carcinoma (HCC)
 Sorafenib response