Standardizing Proteomics Workflow for Liquid Chromatography-Mass Spectrometry: Technical and Statistical Considerations

Introduction: The quantitative measurements based on liquid chromatography (LC) coupled with mass spectrometry (MS) often suffer from the problem of missing values and data heterogeneity from technical variability. We considered a proteomics data set generated from human kidney biopsy material to investigate the technical effects of sample preparation and the quantitative MS. Methods: We studied the effect of tissue storage methods (TSMs) and tissue extraction methods (TEMs) on data analysis. There are two TSMs: frozen (FR) and FFPE (formalin-fixed paraffin embedded); and three TEMs: MAX, TX followed by MAX and SDS followed by MAX. We assessed the impact of different strategies to analyze the data while considering heterogeneity and MVs. We have used analysis of variance (ANOVA) model to study the effects due to various sources of variability. Results and Conclusion: We found that the FFPE TSM is better than the FR TSM. We also found that the one-step TEM (MAX) is better than those of two-steps TEMs. Furthermore, we found the imputation method is a better approach than excluding the proteins with MVs or using unbalanced design.

both unobserved and observed data; (ii) missing at random (MAR) if conditional on the observed data, the MVs are independent of the missing measurements; and (iii) missing not at random (MNAR) when data is neither MCAR nor MAR [30]. The data with missing observations can be analyzed either by excluding the features having missing observations, by using statistical methods that can handle unbalanced data, or by using imputation methods. If the features having missing observations are excluded, then there is loss of information from the experiment. Therefore, the use of methods that can handle MVs, such as imputation methods, are generally preferred. However, the use of imputation methods may lead to wrong interpretation and still these methods are questionable in statistical terms [31,32]. The data set usually consists of biological replicates only or both biological and technical replicates. Biological variability arises from genetic and environmental factors; it is intrinsic to all organisms. The technical approaches include sample collection and storage, sample preparation, extraction, LC separation and MS detection [20]. Sometimes, variations in the biological data or technical approaches to data collection lead to heterogeneity for the samples under study [33,34]. We performed analysis of laser capture microdissection (LCMD)-LCMS high-resolution proteomics dataset using multifactor ANOVA model. We studied the variability in the data based on different TSMs and TEMs. We estimated the contribution of various sources of variation to the overall variability. The study of data variability was done using various analysis methods and transformation and/or normalization techniques. In this paper, we investigated the technical effects of sample preparation and the quantitative MS resulting in heterogeneity for low abundant protein quantification. This will improve the biomarker discovery studies utilizing limited bioreposited tissue resources. We have done all the statistical analysis in R [35] and codes are available from the authors on request.

Proteomics experiment
Data for the methods used in the collection, extraction, and proteomic analysis have previously been published under Hobeika L et al. [36]. Individual data files for MS data (.RAW), peak lists (.mgf), and compressed search results (.mzIdentML) files can be downloaded from the MasslVE data repository (http://massive.ucsd.edu/; MasslVE ID: MSV000079914) and ProteomeXchange data repository [37] (http:// www.proteomexchange.org/; ID: PXD004601). For consideration of variability of the feature detection and MVs abbreviated methods for these studies are provided below.
Tissue collection: FR and FFPE tissue from the same human kidney unsuitable for transplant were cut into 10 μm sections on Polyethylene terephthalate membrane frame slides, stained with Mayer's hematoxylin and glomerular tissue compartments isolated using a Leica LMD6500 Laser Microdissection System.
Protein extraction: Experiments were conducted to compare a single tissue solubilization step using an acid labile surfactant to approaches for tissue decellularization. The single step method used the acid-labile surfactant Protease MAX surfactant with heating (MAX). Two tissue decellularization methods incorporated sequential decellularization with solubilization of the residual pellet with MAX. First tissue decellularization approach used 0.4% SDS + HALT protease/phosphatase inhibitor cocktail (Thermo Fisher) followed by solubilization of residual "ECM" pellet using MAX (SDS.MAX). Second tissue decellularization approach used sequential decellularization with 25mM NH4OH/ 0.5%TritonX-100 (TX) followed by solubilization of residual "ECM" pellet using MAX (TX.MAX). As described in Hobeika L et al. [36], the tryptic peptides were analyzed using a LC-MS Orbitrap ELITE approach with peptide assignments using a Mascot/Sequest search strategy. Scaffold4 was used to set false discovery rate (FDR) control. Finally, we obtained a label-free quantified data of identified proteins (Supplementary File 1). Please see more details about the experimental procedures in "Supplementary File 2". We analyzed the data for comparing statistical methods with MVs in the presence of heterogeneity.

Proteomics data analysis
The purpose of this study is to (1) compare variability between (a) tissue storage methods (TSMs) and (b) tissue extraction methods (TEMs); (2) compare various statistical approaches of analysis and normalization methods.
We have two TSMs (FR and FFPE) and three TEMs (MAX, TX.MAX, SDS.MAX) with three replicates and two MS runs leading to 36 samples (total number of samples = 2 × 3 × 3 × 2 = 36). A flow chart of the experiment is given below in Figure 1.
In the above flowchart, we have shown the basic steps of carrying out the experiment involving TSMs and TEMs. We have repeated the MS two times to get more reliable results for estimating experimental variability. We obtained the following six groups as given below in the Table 1. There are three replicates for each of the six groups thus leading to 18 samples. Then, we have repeated the MS two times for the 18 samples and we obtained six samples for each of the six groups.
Data preprocessing: Initially, there were 728 proteins identified in both runs, 380 proteins identified in run 1 only and 342 proteins identified in run 2 only. There was a total of 1450 identified proteins out of which 1376 proteins were unique, and 37 proteins were redundant and duplicate entries were removed from the data. Furthermore, there were 111 proteins for which all the samples have NA values (MVs). Therefore, we are left with protein data with 1302 proteins that correspond to 1178 gene symbols (Supplementary File 1). The percentage of NA values within each sample (36 samples) ranges from 41.3% −78.3% with a median value of 49.5%. As we have a greater number of groups, therefore it is difficult to perform analysis with this data having MVs. If we discard the proteins having any MVs in any of the samples in a group, then there will be only 26 proteins available. Another way is to retain the proteins having at least one or two observations in each group. A summary of number of proteins available in each group is given below in Table 2. If we use the number of proteins having at least one observation in a group, then we can assess a greater number of proteins. However, we need at least two observations in each group to calculate CV for a protein in each group. Therefore, we used 372 proteins which have at least two observations in each of the six groups for further analysis.
Statistical approaches: The analysis of proteomics data becomes more complex due to non-normality behavior of the data, and greater proportion of MVs within and across the samples. To get a better insight of proteomics data analysis while dealing with these problems, we have performed the analysis using three methods as given below: A1. Method for data excluding missing values: Proteins having complete observations for all the samples, i.e., no MVs, were used for comparison. Proteins having MVs were discarded from the analysis.

A2. Method for data including missing values:
The proteins with MVs across the samples were analyzed using unbalanced ANOVA method [38].
A3. Method for data using imputation: The MVs were imputed after applying the normalization methods to the data [39] as given in next section. We have used the "impute. MAR" function of the R package "imputeLCMD" [40] for imputing the MVs. Three different types of imputation under the assumption of MAR or MCAR, namely, MLE [41], SVD [42] and KNN [43,44] are available in this package. We have used only the SVD method (A3) for imputation.
We applied three different data transformation and/or normalization methods:

N1. Logarithmic transformation:
The raw data is transformed by using logarithmic base 2.

N3. Variance stabilizing normalization:
It is done by applying "justvsn" function available in R package "vsn" [47] to the raw data.
Therefore, by using three methods of analysis (Al, A2 and A3) based on three transformation and/or normalization methods (N1, N2 and N3), we have 9 different combinations (statistical approaches): excluding MVs (A1.N1, A1.N2, A1.N3); including MVs (A2.N1, A2.N2, A2.N3); imputing MVs (A3.N1, A3.N2, A3.N3). We preprocessed the data using these methods to get 9 different datasets (preprocessed data) for 6 groups having 6 samples in each group. We calculated the coefficient of variation (CV) for each protein in  where, y ijk is the transformed and/or normalized data for a protein, α i (i=1,2) is the j th TSM effect,β j , (j = 1,2,3)is the j th TEM effect and (αβ) ij is the interaction effect, TSM×TEM. The term ε ijk is the normally distributed error component and ε ijk N 0, σ 2 . The mapping of the above model to the experimental design allows us to estimate the contribution due to each source of variation for each protein.

Comparison of CV among various groups
We have 141, 372 and 372 proteins obtained by using the analysis methods Al, A2 and A3, respectively. The summary of CV using 9 different statistical approaches for comparisons among TSMs and TEMs is shown below in Table 3. The summary of CV using 9 different statistical approaches for comparisons among six groups of TSM×TEM is shown below in Table 4.

TSM:
We found that median value of CV is lowest in FFPE using all the statistical approaches. Furthermore, within FFPE, the normalization method N3 has the minimum value of median CV for each analysis method. Overall, the minimum median CV is for A1.N3 in FFPE.

TEM:
We have the minimum median value of CV in TX.MAX. We found A1.N2 has the minimum value of median CV.

TSM×TEM:
We have the minimum median value of CV in FR_ MAX followed by FFPE_SDS.MAX using all the approaches. We found A1.N3 has the minimum value of median CV in all the groups except for A1.N2 in FR_SDS.MAX. Overall, the minimum median CV is for A1.N3 in group FR_MAX.
Based on median CV, FFPE is a better choice than FR using all the statistical approaches. Similarly, among TSMs, TX.MAX has the least CV and can be a better choice. However, based on the maximum value of CV, MAX is a better choice for TEM. If we consider approaches (A2 & A3) having greater number of proteins and TEM within FFPE, we see that A3.N3 in FFPE_SDS.MAX is having the least median CV (1.63).

Contribution of Sum of Squares (SS) due to each component
The percent contribution of SS due to each variable to the total SS was computed for each protein. A summary of contribution of each variable to the total variability is given below in Table 5. We found that the TSM has the least contribution to the total variability whereas interaction term has the maximum contribution (SS tsm < SS TEM < SS TSM×TEM ). The imputation method leads to decrease in the SS contribution due to each variable. The proportion of proteins showing significant effects due to TSM, TEM and TSM×TEM using 9 different approaches are given below Table 6. The proportion of proteins showing significant effects due to TSM and TEM and their interaction vary with each statistical approach. The TSM has the least proportion of significant proteins as compared to those of TEM and TSM×TEM. This shows that TSM has the least influence. Furthermore, the imputation approach has the least proportion of significant proteins. This shows that imputation of MVs is a better approach for analysis as it leads to reduction in variability and increase in the number of proteins assessed for analysis.

Analysis for imputed data using VSN
We used ANOVA to test the significance of proteins based on TSM and TEM. The plot of CV of the proteins in increasing order of p-values based on A3.N3 for TSM and TEM are respectively given below in Figures 2 and 3. There are respectively 261 and 296 proteins showing significant effects due to TSM and TEM. From Figure 2, we see that FR has more CV as compared to that of FFPE for most of the proteins. From Figure 3, we found SDS.MAX has more CV as compared to those of MAX and TX.MAX. We applied chisquare test for the proteins having significant effects due to TSM and TEM. We found that there is association between the TSM and the CV (p-value < 0.001). Similarly, in case of TSM, we found that there is association between the variables, TEM and CV (p-value < 0.001). We found that the FFPE is a better method than that of the FR for tissue storage. Further, we found that MAX, the single step approach is better than those of two-step approach for tissue extraction. The maximum contribution to the total variability is due to the interaction effect TSM×TEM and TEM. The TSMs and TEMs have significant effects on the protein expression. However, the effect due to TSM is the least. In the present article, we have used different analysis and normalization methods for the proteomics data. The number of proteins for testing can be increased by either by including the MVs (A2) or by using imputed data (A3). The imputation method (A3) has the least SS contribution than those of A1 (complete data) and A2 (unbalanced data). We found the least proportion of significant proteins when using the imputation method (A3). The normalization method N1, i.e., only logarithmic transformation is not suited for analyzing the proteomics data. The other normalization methods N2 and N3 having lesser CV can be a better approach.

Conclusion
Our study discussed the technical issues with a focus on the statistical analysis. It will provide better insight to the researchers while designing and executing experiments. There may be small changes caused during sample handling and storage, different batches of buffer, electrospray, instrument components, calibration and tuning, etc. While designing any proteomics experiment, we must identify the technical steps with large variability. Therefore, it becomes necessary to understand the data heterogeneity due to biological variability and technical variability of the proteomics methods at each step. We have made the proteomics data available (Supporting file 1). The researchers involved in proteomics research area can use this data for further study. The data can further be used for planning new proteomics experiments. In the future, we will come up with a rigorous statistical approach using different proteomics dataset that could overcome the heterogeneity problem caused due to technical reasons in the proteomics data with MVs. Finally, we can recommend: (i) FFPE is the better choice than FR for tissue storage, (ii) one-step TEM is better than the two-step TEM, (iii) Imputation method (A3) is the best approach, (iv) N2 or N3 method of normalization should be the preferred choice.

Supplementary Material
Refer to Web version on PubMed Central for supplementary material.

Acknowledgement
Authors wish to acknowledge the Indian Council of Agricultural Research, India and the University of Louisville, United States of America for providing facilities to carry out the research.

Funding
This work was partially supported by the National Institute of Environmental Health Sciences (P42ES023716-PI: Dr. Sanjay Srivastava) and the National Institute of Health (P20 GM113226-PI: Dr. Craig McClain).

List of Abbreviations:
LC Liquid Chromatography    Summary of number of proteins and missing values in different groups.   Summary of CV using 9 statistical approaches among six groups of TSM×TEM.