PMCCPMCCPMCC

Search tips
Search criteria 

Advanced


 
Formats
Article sections
Authors
Related links
Logo of protsciBioMed CentralBiomed Central Web Sitesearchsubmit a manuscriptregisterthis articleProteome ScienceJournal Front Page
 
Proteome Sci. 2011; 9: 73.
Published online 2011 November 17. doi:  10.1186/1477-5956-9-73
PMCID: PMC3247853
Copyright ©2011 Kato et al; licensee BioMed Central Ltd.
Variance decomposition of protein profiles from antibody arrays using a longitudinal twin model
Bernet S Kato,corresponding author1,5 George Nicholson,2 Maja Neiman,3 Mattias Rantalainen,2 Chris C Holmes,2 Amy Barrett,4 Mathias Uhlén,3 Peter Nilsson,3 Tim D Spector,1 and Jochen M Schwenkcorresponding author3
1Department of Twin Research and Genetic Epidemiology, King's College London, St. Thomas' Hospital, Westminster Bridge Road, London SE1 7EH, UK
2Department of Statistics, University of Oxford, 1 South Parks Road, Oxford OX1 3TG, UK
3Science for Life Laboratory Stockholm, KTH - Royal Institute of Technology, Box 1031, SE-171 21 Solna, Sweden
4Oxford Centre for Diabetes, Endocrinology and Metabolism, Churchill Hospital, Old Road, Headington, Oxford, OX3 7LJ, UK
5Respiratory Epidemiology and Public Health, Imperial College London, Manresa Road, London SW3 6LR, UK
corresponding authorCorresponding author.
Bernet S Kato: b.kato/at/imperial.ac.uk; George Nicholson: nicholso/at/stats.ox.ac.uk; Maja Neiman: maja.neiman/at/scilifelab.se; Mattias Rantalainen: rantalai/at/stats.ox.ac.uk; Chris C Holmes: cholmes/at/stats.ox.ac.uk; Amy Barrett: amy.barrett/at/drl.ox.ac.uk; Mathias Uhlén: mathias.uhlen/at/scilifelab.se; Peter Nilsson: peter.nilsson/at/scilifelab.se; Tim D Spector: tim.spector/at/kcl.ac.uk; Jochen M Schwenk: jochen.schwenk/at/scilifelab.se
Received August 8, 2011; Accepted November 17, 2011.
This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
Abstract
Background
The advent of affinity-based proteomics technologies for global protein profiling provides the prospect of finding new molecular biomarkers for common, multifactorial disorders. The molecular phenotypes obtained from studies on such platforms are driven by multiple sources, including genetic, environmental, and experimental components. In characterizing the contribution of different sources of variation to the measured phenotypes, the aim is to facilitate the design and interpretation of future biomedical studies employing exploratory and multiplexed technologies. Thus, biometrical genetic modelling of twin or other family data can be used to decompose the variation underlying a phenotype into biological and experimental components.
Results
Using antibody suspension bead arrays and antibodies from the Human Protein Atlas, we study unfractionated serum from a longitudinal study on 154 twins. In this study, we provide a detailed description of how the variation in a molecular phenotype in terms of protein profile can be decomposed into familial i.e. genetic and common environmental; individual environmental, short-term biological and experimental components. The results show that across 69 antibodies analyzed in the study, the median proportion of the total variation explained by familial sources is 12% (IQR 1-22%), and the median proportion of the total variation attributable to experimental sources is 63% (IQR 53-72%).
Conclusion
The variability analysis of antibody arrays highlights the importance to consider variability components and their relative contributions when designing and evaluating studies for biomarker discoveries with exploratory, high-throughput and multiplexed methods.
Keywords: Variance decomposition, linear mixed-effects model, longitudinal twin study, suspension bead arrays, antibodies, protein profiling
Introduction
There is an enormous unmet need for biomarkers to characterize disease type, status, progression, and response to therapy. A biomarker, which can be defined as a physical sign or laboratory measurement that serves as an indicator for biological processes, pathogenic processes, or pharmacologic responses to a therapeutic intervention [1,2], can then potentially be used for screening, prognosis, monitoring response to treatment, and detection of recurrent diseases. Technological advances in the different fields of life science have enabled the generation of data from high-throughput screening experiments, which can then be used to identify novel biomarker molecules.
Proteomics, the large-scale study of the protein content of cells, organs or organisms, offers the potential to evaluate global changes in protein expression, protein interaction patterns and their post-translational modifications that take place in response to normal or pathological stimuli. The availability of DNA microarray analysis permits expression of thousands of genes to be monitored simultaneously, but RNA expression has been found to correlate significantly to only one third with the corresponding actual protein content [3]. Therefore the importance of proteomic research cannot be overstated, as proteins within the cell provide a structural and functional framework by producing energy and allowing communication, movement, and reproduction [4].
Recent technological advances in the field of genomics and proteomics have also brought about a new statistical research area, the analysis of data from high-throughput screening experiments. Using protein microarray technology, the profiles of thousands of proteins can be monitored in a high-throughput manner with the aim of selecting a subset of proteins that could be characterised as potential biomarkers creating possibilities for diagnostic, prognostic and disease progression monitoring [5]. When proteomics emerged, scientists hoped that available technologies would help them sift through the proteome of blood and tissues to identify profiles of proteins that indicated the presence of a disease. Researchers could then use these molecular signatures as biomarkers to diagnose diseases in their formative stages and therefore be able to tailor appropriate intervention and treatment. However, one problem is that the most preferred sample types, serum or plasma, are highly complex and influenced by a protein composition where 99% of the total protein mass is covered by only 20 proteins. Even though all proteins are likely to appear in blood at a given point in time, due to function, secretion or leakage, their levels can vary by more than 10 orders of magnitude, with abundant proteins masking the presence of rare ones. Among the proteomic approaches applied in the quest for novel biomarkers today, mass spectrometry, often in combination with 2D-gel electrophoresis or chromatographic separation, is still the most commonly used technique in the field [6], even though the technique is limited with regard to new diagnostics [7]. Affinity-based alternatives are on the rise to enable proteome-wide analysis, and initiatives have been set up to produce, validate and offer a larger number of validated affinity reagents [8]. Platforms such as microarrays provide solutions for the implementation of larger sets of affinity reagents and analysis of multiple parameters on a minute amount of sample within a single experiment [9]. Among the immunoassay platforms that employ directly labelled samples [10-12], antibody suspension bead arrays to perform highly multiplexed protein profiling in a large number of serum samples [13].
Utility of twin studies
The phenotype of an individual might change as a result of complex, dynamic interactions between an individual's genome and environmental factors. To characterize and quantify the contributions that genes, the shared and the individual-specific environment and the interactions of these make to human complex traits and phenotypes, we use data from relatives who grow up in similar environments but are of differing genetic relatedness (the so-called 'twin design'). For a long time, twin studies [14] have been a valuable source of information about the genetic basis of complex traits. Twins provide a powerful design for inferring genetic effects, as they are blocked for in utero, dietary and socio-economic effects due to common upbringing. Identical/monozygotic (MZ) twins derive from a single fertilized egg and therefore inherit identical genetic material while non identical/dizygotic (DZ) twins share 50% of their genetic material. In a classical twin study the assumption is that MZ and DZ twins share the same degree of similarity in their environments and therefore one compares phenotypic concordance within MZ twin pairs to the concordance within DZ twin pairs. Comparing the resemblance of MZ twins for a trait or disease with the resemblance of DZ twins offers the first estimate of the extent to which genetic variation determines phenotypic variation of the trait or disease. The sources of variation commonly considered in genetics are additive genetic influences (sum of individual effects of all loci that influence the phenotype), dominant genetic influences (these represent interactions between alleles at all loci that influence the phenotype), common environment (influences shared by family members e.g. socio-economic status) and individual environment (influences that cause differences among members of one family, e.g. dietary and lifestyle choices that are not shared by other family members. Two important population parameters are heritability and familiality. The narrow-sense heritability is the proportion of phenotypic variation within a population attributable to additive genetic effects. Familiality is the proportion of phenotypic variation within a population attributable to genetic (additive and dominant) and common environmental effects.
In the following, we have investigated the variance of exploratory affinity arrays using longitudinal twin study.
Results
Important attributes of discovered biomarkers, in addition to being associated with a disease, are that they can be measured with precision (low variance) and exhibit relatively low amounts of short-term variability relative to the effect size of the difference between healthy and diseased individuals [15]. For antibody-based proteomics, antibodies that have high proportions of variability attributable to familiality and individual environment and have relatively small technical variability are therefore potentially more valuable for use in biomarker discovery.
Pre-processing
For this study the measured relative quantity of interest is the intensity of fluorescence emitted and measured when antibodies immobilized onto beads capture their target protein from a sample (see Figure ​Figure1).1). The intensity levels are influenced by a number of underlying factors besides the antigen's abundance. These are antibody-antigen kinetics and affinities, the number of accessible binding sites on both immobilized antibodies and the labelled antigens, the number of biotin molecules per antigen being incorporated during labelling reaction, as well as complex formation introducing additional biotin molecules carried by proteins attached to the antigen of interest. Therefore the variation in the measured intensity levels will be influenced by both biological and experimental factors including sample preparation, incubation and read-out.
Figure 1
Figure 1
Experimental workflow - The process begins with (1.) the distribution of samples into microtiter plates according to a defined layout, dilution and heat treatment. (2.) The protein content of the samples is then label with biotin and (3.) the samples (more ...)
Molecular phenotypes can be noisy, context dependent, and sensitive to laboratory, equipment and experimental conditions. On the current platform, the fluorescence intensities obtained on each sample can in general not be compared directly due to variations in sample treatment, labelling and detection, thus appropriate pre-processing is needed. Plots of the signal intensities of each antibody across the randomized samples (Figure ​(Figure2A2A and Additional File 1) show that for some of the antibodies, such as for example apoh-HPA001654, renbp-HPA000428, apoh-HPA003732 and icam1-HPA004877, signal intensities tended to decrease as intensities were measured from plate well 1 through plate well 96 suggesting an 'intensity-drift effect'. This could be introduced during sample preparation, modification of the potential binding site during sample labelling, slight bleaching of the fluorophore of the reporter protein, or by the dissociation of antibody-antigen complexes over time the samples are being measured.
Figure 2
Figure 2
Protein profiles. A) Profiles before any data processing are shown from four antibodies across all samples. B) Profiles for the antibodies are shown after normalization, removal of outliers and BoxCox transformation. The red line indicates the locally (more ...)
The data were normalized using probabilistic quotient normalization (PQN)[16]. To assess the efficacy of the normalization we looked at the within-sample concordance before and after normalization, by computing the Spearman correlation between the 48 pairs of aliquot replicates for each antibody in the data set. See appendix for details on the correlation model and an example in Figure ​Figure3.3. Across the 69 antibodies, before normalization the median Spearman correlation was 0.26 with interquartile range (IQR 0.19 - 0.36) and after normalization, median Spearman correlation was 0.36 (IQR 0.24 - 0.44). These results show that normalization improves the within-sample concordance.
Figure 3
Figure 3
Correlation model. Suppose we have three duplicated samples 1, 2 and 3 assayed at 6 antibodies B1 - B6. Denoting the aliquot pairs of the samples as 1a, 1b, 2a, 2b, 3a and 3c (see table above), to check the within sample concordance we look at correlations (more ...)
Outliers
Using singular value decomposition (SVD) [17] we projected the 69 dimensional antibody space into a two-dimensional space (two principal components). Prior to performing SVD on the data matrix, we scaled the data of each variable (antibody) to have a mean of 0 and standard deviation 1. Plots of the two principal components are shown in Figure ​Figure4.4. As can be seen in the figure, most of the samples cluster together. However, there are some outlying samples. For each of the sample points in the figures, Mahalanobis [18] distance from the origin was calculated. Assuming that the principal components have a bivariate normal distribution, the Mahalanobis distance of each sample was compared to a chi-squared distribution with 2 degrees of freedom. A sample was defined as an outlier if it had a Mahalanobis distance bigger than 18.42 (critical value at a 0.0001 level of significance). Based on the above criteria, the data set contains 6 out of 270 samples classified as outliers as shown in Figure ​Figure44 highlighted using numbers.
Figure 4
Figure 4
Principal component projection. The graph shows the normalized data for the 270 samples in the dataset in the two-dimensional antibodies space. Outlying samples are marked with numbers.
Variance decomposition
Plots of the data after normalization, removal of outliers and Box-Cox transformations are presented in Figure ​Figure2B2B and Additional File 2. Further, table ​table11 shows for each antibody the AIC obtained from fitting models 1 and 2 (see Modelling and Estimation in the Materials and Methods Section) on these data. As can be seen in the table, for some antibodies model 1 had a better fit (smaller AIC), whilst for others model 2 had a better fit. This shows that for some antibodies including a linear plate specific intensity-drift effect modelled the data better. Subsequently for each antibody, model 1 or 2 (depending on the model that had a smaller AIC) was used to obtain maximum likelihood estimates of the fixed effects and the variance components.
Table 1
Table 1
Akaike's information criterion
For each antibody the percentage of total variance attributable to familiality (fam), individual environment (env), common visit (cv), individual visit (iv) and residual experimental effects (exp) each antibody is presented in Table ​Table22 (rows are ordered by decreasing familiality). In Table ​Table2,2, we observe that for many proteins the variances of the factors other than experimental variance are estimated to be zero. In general, the maximum likelihood (ML) estimators of the variance components often lie on the boundary of the parameter space i.e. close to zero. This can result from (a) the underlying "true" values of the variances being relatively close to zero and/or (b) there being less information about these variance components, implying that the ML estimators have relatively high variance.
Table 2
Table 2
Phenotypic variance
Figure ​Figure5A5A represents bar plots of the decomposition of the total phenotypic variance for each antibody after normalization, removal of outliers and Box-Cox transformations. The graph in Figure ​Figure5B5B summarizes the decomposition of total phenotypic variance across all antibodies. As can be seen from these figures, most of the antibody-derived intensities have a high residual (experimental) variance suggesting that a lot of the variance is attributable to experimental sources and sample complexity. Across the 69 antibodies, the median proportion of total phenotypic variance attributable to familiality is 12 (IQR 1-22%); the median proportion of total phenotypic variance attributable to experimental sources is 63% (IQR 53-72%). Combining familiality and individual environment renders the proportion of non-experimental variance that is stable over the time period between the two visits during which data was collected. Across, the 69 antibodies, the median proportion of total variance attributable to familiality and individual environment is 25% (IQR 14-33%). The median proportions of total phenotypic variance attributable to individual visit and common visit is 6% (IQR 0-18%) and 0% (IQR 0-3%), respectively.
Figure 5
Figure 5
Variance decomposition. A) The bar plot summarizes the decomposition of total phenotypic variance of each antibody. The colours in each bar represent the proportion of total phenotypic variability attributable to familiality (fam), individual environment (more ...)
In what follows, we present summaries of the proportion of non-experimental variance attributable to familiality, individual environment, common visit and individual visit. The boxplots in Figure ​Figure5C5C summarize the decomposition of non-experimental variance across all antibodies. Across the 69 antibodies the median proportion of non-experimental variance attributable to familial sources is 34% with (IQR 3-60%).. Across all 69 antibodies, the median proportion of non-experimental variance attributable to familiality and individual environment is 71% (IQR 51-93%). The median proportions of non-experimental variance attributable to common visit and individual visit are 0% (IQR 0-10%) and 16% (IQR 0-40%), respectively, indicating that most of the unstable (short-term dynamic) biological variation is attributable to individual visit effects.
Antibodies for targeting potential biomarkers
Highlighting those HPA antibodies with most preferred characteristics, we identified five interesting antibodies where more than 50% of the total phenotypic variance is attributable to stable biological (familial and individual environmental) variability. Among these target proteins were interleukin 12 (il12a), a heterodimeric 70 kDa glycoprotein with implications in cell-mediated immunity [19] targeted by HPA001886 (62%), and angiotensin-converting enzyme 2 (ace2) targeted by HPA000288 (59%). Ace2 is an exopeptidase that catalyses the conversion of angiotensin I or II and its pharmacological inhibition is associated with protective effects from cardiovascular diseases. Further, ADAMTS-like protein 4 (adamstsl4) targeted by HPA006279 (57%) is a secreted glycoprotein with a suggested function to be involved in cell adhesion and protease regulation [20]. The tissue-type plasminogen activator (plat) targeted by HPA003412 (54%) is a secreted serine protease which converts the proenzyme plasminogen to plasmin and up-regulated gene expression has been associated with heart failure [21]. The transforming growth factor beta receptor III (tgfbr3) targeted by HPA008257 (53%) is a glycosylated protein that is found as membrane-bound and as soluble form, and is among the most widely distributed TGF-beta family members and is suggested to function as accessory molecules in TGF-beta and FGF systems [22]. Interestingly, a strong negative correlation between profiles of interleukin 12 and tissue-type plasminogen activator was observed (Figure ​(Figure6).6). A linkage between these two proteins has been reported before [23] with il12a inducing the activation of fibrinolysis and thrombin generation. It was also shown that plasma levels of plat decrease shortly after il12 administration but increased over time as il12a levels lowered.
Figure 6
Figure 6
Correlation of profiles of HPA001886 and HPA003412. A) The intensity profiles from antibodies targeting IL-12 and PLAT were found to be strongly negatively correlated (R = -0.95, red line), which suggested a biological connection between these two proteins. (more ...)
The results suggest that some antibodies are more suitable and promising than others and that they can be potential candidates to be used further in biomarker discovery studies. Among these, HPA006279 has been shown to reveal interesting differences in plasma protein profiles in Metabolic syndrome cases and controls [24]. The identified antibodies should now be further investigated (i) separately, in combination or supplemented by other candidates in high-throughput screenings across diseases, and (ii) to verify their potential implications in diseases related cardiovascular disorders by meeting the criteria of clinical test systems (sandwich assays), which will most likely the improve precision (reduce variance).
Discussion
In this paper we use a linear mixed-effects model in a unique twin design with duplicated repeated measurements to apportion the total variation of molecular phenotypes (protein profiles) into biological and experimental (technical) variation. Understanding the sources of variation (e.g familial, individual environmental, and experimental) inherent in the measurement of a molecular phenotype is a key step in assessing the potential for stable, informative biomarkers. We observed that across the 69 antibodies the median proportion of total variation attributable to familial sources was 12%. This familial component is consistent with protein profiles having the potential to reflect the polygenic basis of complex disease. Further, across the 66 HPA antibodies, the median proportion of total variation explained by stable sources (familiality and individual environment) was 25%. Stable variation in the current setting comprises that a protein profile remains constant within an individual over the course of the sampling period. A small proportion of variation originated from the short-term biological component, individual visit (median proportion 6%). Common visit represented an inconsiderable amount of variation. Most of the variation originated from experimental sources (median proportion 63%). To our current knowledge, the present study is the first to address the key issue of investigating sources of variation in data generated by exploratory antibody microarrays. It should provide important information when aiming at designing and utilizing such assays and be valuable for multiplexed and quality controlled assays [25] that will become more widely used and accepted for clinical testing.
Ultimately, diagnostic tools built on markers discovered via these screening approaches could become valuable approaches to predict disease state and progression. As an example, antibodies that allow the comparison of individuals that are discordant for diabetes and diabetes-related clinical traits would be useful for identifying individuals likely to suffer from diabetes in the future, long before conventional diagnostic techniques can prove effective. In this way these affinity-based proteomics discoveries would become useful in clinical settings.
However, most antibodies used in this screening method have a large residual variance suggesting that a large proportion of variation in the data is experimentally derived. Potential (and non-separable) sources of this experimental variation, which exclude sample collection, preparation and storage, are the:
(i) complexity and composition of the serum samples which has an effect on the assays;
(ii) biotin modification of samples with regard to the numbers and variability of modifications introduced per molecule and sample;
(iii) sample treatment in terms of liquid transfers, heating, and assay buffer dilution;
(iv) assay procedure with immobilized antibodies selectively capturing aggregated or free molecules from the surrounding solution;
(v) fluorescence-based read-out being influenced by bleaching and dye incorporation onto the target molecules.
(vi) specificity of antibody binding events.
Addressing these issues would lead to a reduction in the proportion of phenotypic variation arising from experimental sources. This would in turn reduce the sample sizes (or degree of technical replication) required to detect epidemiological effects of interest. A reduction of the experimental variability, possible to achieve by using two antibodies to detect a target protein, would ensure that the experimental noise does not swamp biological signals of interest. The technical precision of such a measurement and of other antibodies can be improved, either by addressing the issues outlined above, or, in the short term, by assaying samples in technical replicate.
Research in the field of proteomics is advancing, with affinity-based approaches emerging alongside classical mass spectrometric approaches. With array-based proteomics becoming a promising area in the field of biomedical research, decomposition of the underlying variation in protein profiles into biological (both stable and longitudinally fluctuating) and experimental components is an important and useful step in exploring the applicability of antibody arrays for the exploration of the proteome. Ultimately, such proteomic strategies may lead to new disease markers and drug targets can be identified, benefitting from the possibilities offered by the versatility of both the employed affinity reagents and multiplexed techniques.
Materials and methods
In this paper we longitudinally collect serum samples from a cohort of 77 twins (56 MZ pairs and 21 DZ pairs) to explore the sources of variation underlying antibody array-based protein phenotypes with the aim of aiding the design and interpretation of future proteomic research. We use a mixed-effects modelling approach to obtain estimates of variance components.
This is a longitudinal family study whose design involves duplicate measurements over time in a sample that includes related individuals (twins). Consequently it combines the features of a longitudinal study (the power to study the pattern of changes in a trait over time) with the features of a cross-sectional twin study (which allows one to estimate the variation in a trait attributable to familial sources).). The study design allows us to break down the observed molecular phenotypic variance into four biological components. These are familiality (the combined effects of genetic and common environment) and individual environment, both of which are longitudinally stable biological components. The two short-term biological components are common visit and individual visit, which respectively measure the amount of shared (by twins within a pair) and non-shared longitudinal variation observed over the sampling period of the present study (i.e several months).
Subject recruitment
A total of 154 healthy, postmenopausal, female twins, comprising of 56 MZ and 21 DZ pairs, were ascertained from the TwinsUK database at St. Thomas' Hospital http://www.twinsUK.ac.uk, and invited to participate in this study. Recruitment of the twins was as follows: eligible twins were sent an information sheet containing details of the study as well as two consent forms. If the twins agreed to participate they signed and sent one copy of their consent forms to the administration team at the Department of Twin Research and genetic epidemiology (DTR), King's College London. Once a consent form was received, a member of the administration team would contact the twins by letter and phone to book their appointment. In addition, 68 twins (34 MZ pairs) were asked to re-attend once they had successfully completed the first visit with the time between the two visits ranging from 63 to 99 days (26 twins), 104 to 140 days (34 twins) and 150 to 216 days (8 twins). At each visit to the clinic, the twins were in pairs.
Fasting blood samples were collected from all selected twins from one or both visits rendering a total of 222 samples with samples from individuals of a twin pair taken on the same day. These were to be used for the study of biological and technical variability in the various molecular phenotyping platforms, focusing on the unique opportunities afforded by studies in a large set of clinically well characterised monozygotic and dizygotic twins. These studies were all part of the Molecular Phenotyping to Accelerate Genomic Epidemiology (MolPAGE) project. For the MolPAGE project, eligible volunteers were healthy, Caucasian females of North Europe descent, aged between 45-76 years old. The study was approved by St. Thomas' Hospital Research Ethics Committee (EC04/015 Twins UK). The experiments and subsequent analyses presented in this paper are based on serum samples extracted from the fasting blood samples.
Antibodies, Bead Coupling, Sample preparation and Assay procedure
Antibodies
Protocols for antigen selection, cloning, protein expression, immunization of rabbits, and affinity purification to yield polyclonal antibodies were performed as described previously [26,27]. All protein fragments used for immunization were produced with a tag and a target protein part of on average 120 amino acids and in general those proteins had a size of approximately 30 kDa. Antibodies from the Human Protein Atlas (HPA) have been selected as described elsewhere [24] (see also Table ​Table1)1) and 66 were involved the study, all tested on protein arrays for specificity [28]. In addition rabbit IgG (Jackson ImmunoResearch), a human serum albumin binding protein (Affibody AB) and an anti-IgA antibody (Dako Cytomation) were also included but not highlighted in further detail.
Bead coupling
Antibodies were coupled to magnetic carboxylated beads (MagPlex Microspheres, Luminex-Corp.) in accordance to the manufacturer's protocol with minor modifications. For each antibody, 3.2 μg were coupled to 106 beads in a 96-well half-area microtiter plate (Greiner BioOne) and beads were allowed to sediment for 30 s on a magnetic bead separator (LifeSep™, Dexter Magnetic Technologies Inc.) prior to the removal of supernatants. For an albumin binding Affibody® molecule, 4 μg were coupled. Beads were stored in a protein-containing buffer (Blocking Reagent for ELISA, Roche) with NaN3. All coupled beads were re-suspended with sonication in an ultrasonic cleaner (Branson Ultrasonic Corporation) for 5 min prior to storage at 4°C. All antibody-coupled beads were counted using the Luminex LX200 instrument and the coupling efficiency for each antibody was determined via R-Phycoerythrin labeled anti-rabbit IgG antibody (Jackson ImmunoResearch). A bead mixture was created in solution, and optimized as previously described [29].
Sample labelling
Protocols for sample labelling and assay procedure were based on earlier work [13]. At first, samples were thawed at room temperature and centrifuged for 10 min at 10,000 rpm, transferred into a microtiter plate (Abgene) and heat treated at 56°C over 30 min in a thermo cycler (Biorad). The plates were centrifuged (1 min at 3,000 rpm) and 3 μl of each sample was added to 24.5 μl 1 × PBS with a liquid handler (Plate mate 2 × 2, Matrix). N-hydroxysuccinimidyl ester of biotinoyl-tetraoxapentadecanoic acid (NHS-PEO4-Biotin, Pierce) was then added at 10-fold molar excess to yield an overall 1/10 sample dilution followed by an incubation over 2 hours at 4°C in a microtiter plate shaker (Thermomixer, Eppendorf). The reaction was stopped by the addition of Tris-HCl, pH 8.0 at a 250-fold molar excess over biotin and incubated for another 20 min at 4°C prior to a final storage at -80°C.
Assay procedure
All samples were utilized without removing unincorporated biotin and diluted 1/50 utilizing a liquid handler in a assay buffer composed of 0.5% (w/v) polyvinylalcohol and 0.8% (w/v) polyvinylpyrrolidone (Sigma) in 0.1% casein in PBS (PVXC) supplemented with 0.5 mg/ml non-specific rabbit IgG (Bethyl). A second heat treatment was performed as above and 45 μl sample were added to 5 μl of bead mixtures in a microtiter plate (Greiner BioOne). Incubation took place over night on a shaker at an ambient temperature. Beads were washed in wells with 3 × 50 μl PBST (1 × PBS pH 7.4, 0.1% Tween20) on a magnetic bead separator with a liquid handling system followed 10 min with 50 μl of a stop solution containing 0.1% paraformaldehyde in PBS. Beads were washed 1 × 50 μl PBST and 30 μl of 0.5 μg/ml R-Phycoerythrin labeled streptavidin (Invitrogen) in PBST were added and incubated for 20 min. Finally, beads were washed 3 × 50 μl and measured in 100 μl PBST. Measurements were performed using a Luminex LX200 instrument with Luminex xPONENT software and counting at least 50 events per bead ID, displaying binding events as median fluorescence intensity. A summary of the experimental work flow is presented in Figure ​Figure1.1. In what follows we will refer to each measured proteomic phenotype in the study as an "antibody".
Of the 222 serum samples obtained from fasting blood, 48 samples from 24 MZ twin pairs were split into two aliquots (which usefully allows assessment of technical replicability) rendering an additional 48 sample aliquots. The resulting 270 sample aliquots were then randomized and placed onto three 96-well plates, where each plate contained an extra 6 reference samples. Subsequently a total of 69 antibodies were utilized to render a 270 by 69 matrix of protein (antibody) array data.
Statistical analysis - Data pre-processing, Modelling and Estimation
Data pre-processing
Probabilistic quotient normalization (PQN) [16] was used to pre-process the data. In metabonomics, PQN has been used as a robust method to account for dilution of complex biological mixtures. The method performs well in compensating for different dilutions in samples. We found this normalization method to increase a correlation-based measure of technical reproducibility across two aliquots of the same sample (see Results).
Outliers
Intensity distributions of the antibodies after normalization using PQN were investigated to identify samples that were potential outliers. In order to detect potential outliers we used singular value decomposition [17] to project the data into two-dimensional space. Projecting the data onto the space spanned by the first two principal-component loadings vectors is a commonly used technique for detecting anomalous observations whose extreme behaviour explains a considerable proportion of variance in the data.
Modelling and Estimation
All subsequent analyses were done on data that had been normalized using PQN and had undergone removal of outlying samples. Further, a Box-Cox transformation [17] was performed on the data from each antibody. This transformed the data so that their distribution more closely resembled a Gaussian, thereby increasing the quality of fit to the data of the Gaussian-based mixed-effects model.
Variability analysis
For each molecular phenotype (antibody) we used a statistical model to quantify the biological and experimental variation inherent in the phenotype. As mentioned before, the aim was to partition the total variability in a molecular phenotype into variability that is attributable to familiality (additive genetic, dominant genetic and common environmental), individual environmental, common visit, individual visit and experimental components, respectively. Note that familial and individual environmental variability together make up stable biological variability.
As mentioned previously, the measured quantity of interest was the intensity of fluorescence emitted when an immobilised antibody captured an antigen in a sample. For each antibody let Yijkl denote the normalized and transformed intensity of the lth aliquot replicate (1, 2), at the kth visit (1, 2) of the jth twin (1, 2) from the ith twin pair (1, 2, ...,77). The data for each antibody was modelled using a linear mixed-effects model [30-32]:
equation M1
(1)
where μp(i, j, k, l) is the mean intensity for samples on plate p (the function p(i, j, k, l) maps sample aliquots to plates), A is the additive genetic effect, D is the dominant genetic effect, C is the common environment effect, E is the individual or unique environment effect, W is the common visit effect, V is the individual visit effect, and lastly ε represents the residual experimental effects. The plate mean μp(i, j, k, l) is a fixed effect, while A, D, C, E, W, and V are random effects. It is assumed that the random effects are additive and independent. This implies that the total phenotypic variance can be decomposed as:
equation M2
The covariance between phenotypic measurements on a pair of samples gathered at the same visit from monozygotic twins is Var(A) + Var(D) + Var(C) + Var(W), and for dizygotic twins is (1/2)Var(A) + (1/4)Var(D) + Var(C) + Var(W). Familiality (fam) is the proportion of variance attributable to genetic and common environmental effects:
equation M3
For twin data, the variances Var(A), Var(C) and Var(D) are not all identifiable. However familiality is an estimable function of these variance parameters.
For each antibody, maximum likelihood (ML) estimates of the variance components were obtained by fitting model (1) using the lmer function from the lme4 package in R [33] after normalization, removal of outliers and Box-Cox transformations. The model was reduced prior to fitting. The variances of A, C and D are not identifiable. In order to address non-identifiability we re-parameterized the three non-identifiable variances, Var(A), Var(C) and Var(D) by defining two new random effects H and M (whose variances are identifiable) as follows:
equation M4
equation M5
Subsequently the familial variance (which represents the combined effects of genetics and common environment) was obtained as the sum of the estimated variances of H and M, i.e familial variance = Var(H) + Var(M).
Note that if x1 and x2 are the phenotype measurements of two monozygotic (MZ) twins then the covariance matrix of the measurements is such that: Var(x1) = Var(x2) = Var(H) + var(M) +Var(E) and covariance(x1, x2) = Var(H) + Var(M). Similarly if x1 and x2 are the phenotype measurements of two dizygotic (DZ) twins then covariance matrix of the measurements is such that: Var(x1) = Var(x2) = Var(H) + Var(M) + Var(E) and covariance(x1, x2) = Var(H). A similar approach has previously been used and elaborated on metabolomics data [34].
The PQN pre-processing procedure that we applied should mitigate any systematic differences between measurements in the different wells introduced by the instrumental time drift between plate well 1 through to plate well 96. However, in order to investigate if there were any substantial significant intensity-drift effects remaining after PQN we also fit a second model (hereafter referred to as model 2) by introducing a linear, plate-specific well location term into model 1, where well location is a number between 1 and 96 representing the position of a sample in the plate wells. That is in model 2, for each antibody we fit a linear mixed-effects model with linear, plate-specific intensity-drift fixed effects and the random effects A, D, C, E, W, and V as before. Subsequently Akaike's information criterion [35,36] was used to compare the fit of models 1 and 2. Akaike's information criterion (AIC) is a measure of the goodness of fit of an estimated statistical model and is defined as:
equation M6
where k is the number of parameters in the model, and L is the maximized value of the likelihood function for the estimated model. For a given a data set, several competing models may be ranked according to their AICs and the model with the smallest AIC is chosen as the one that fits the data best. The pre-processed data and R script used to perform the variability analysis can be found as Additional Files 3 and 4.
Abbreviations
HPA: Human Protein Atlas; IQR: Interquartile range; Fam: familial sources; Env: individual environment effects; Cv: common visit effects; Iv: individual visit effects; Exp: experimental effects; MZ: monozygotic; DZ: dizygotic; PQN: probabilistic quotient normalization; AIC: Akaike's information criterion; LOWESS: locally weighted scatterplot smoothing.
Competing interests
The authors declare that they have no competing interests.
Authors' contributions
BSK did the data analysis and generated the manuscript, GN provided guidance on the data analysis and helped generate the manuscript, MN has carried out array experiments and participated in data analysis, MR provided guidance on the data pre-processing and analysis, AB was responsible for storage and distribution of the serum samples, CCH supervised aspects of this manuscript, MU supervised aspects of this manuscript, PN has edited the manuscript and supervised the experiments, TDS supervised aspects of this manuscript, JMS generated the manuscript, designed and carried out array experiments and participated in data analysis. All authors read and approved the final manuscript.
Appendix
Correlation model
Consider a proteomics platform, which measures the abundance of p proteins, indexed by k [set membership] {1,...,p}. Suppose that there are n samples, indexed by i [set membership] {1,...,n}, and each is measured in technical duplicate, with duplicates indexed by j [set membership] {1, 2}. A model for the observed data is:
equation M7
where:
• equation M8is the abundance of the kth protein as measured in the jth replicate of the ith biological sample
• μ(k) is the mean (across all samples) of the measured abundance of protein k
• equation M9 is the biological deviation in the abundance of protein k in sample i around μ(k)
• equation M10 is the experimental variation in replicate j of sample i at protein k.
The correlation-based measure of technical reproducibility used in the current paper was calculated for each protein in turn, as follows. For the kth protein, we formed the two vectors equation M11 and equation M12, and then estimated the correlation between the two: equation M13; r(k) is directly related to the utility of the platform, as it measures the reproducibility across replicates for each individual protein.
A different, commonly used (and often misinterpreted) measure of correlation would be calculated as follows. For the pair of replicates of sample 1, form the two vectors equation M14 and equation M15, and then estimate the correlation between the two: s1 [equivalent] cor(w11,w12). It is less appropriate to use the correlation across proteins, s1, as this is affected by the overall scale of the range of measurement. For instance, for the same level of technical variation, equation M16, s1 can be made arbitrarily large by increasing the range between the expected values of the lowest- and highest-expressed protein levels (i.e. increasing the range of the μ(k)).
Supplementary Material
Additional file 1
Protein profiles before data processing. Profiles from all antibodies across all samples are shown before any data processing, with the red line indicating the locally weighted scatterplot smoothing (LOWESS).
Click here for file(272K, PDF)
Additional file 2
Protein profiles after data processing. Profiles from all antibodies across all samples are shown after any data processing, with the red line indicating the locally weighted scatterplot smoothing (LOWESS).
Click here for file(283K, PDF)
Additional file 3
Pre-processed data. The csv file contains pre-processed intensity values for all antibodies across all samples.
Click here for file(169K, CSV)
Additional file 4
R script for variability analysis. The R script can be used to perform the variability analysis with the pre-processed data from Additional File 3.
Click here for file(5.1K, R)
Acknowledgements
This work is part of MolPAGE, the Molecular Phenotyping to Accelerate Genomic Epidemiology project (European Union grant LSHG-512066, 6th framework funding programme). We would gratefully like to thank Mark I McCarthy, John Bell, and Maxine Allen for coordinating this project. For the valuable comments on this manuscript we like to thank Krina Zondervan and Kourosh Ahmadi. GN and CCH acknowledge funding from MRC Harwell, UK. CCH acknowledges funding from the Wellcome Trust.
References
  • Lesko LJ, Atkinson AJ Jr. Use of biomarkers and surrogate endpoints in drug development and regulatory decision making: criteria, validation, strategies. Annu Rev Pharmacol Toxicol. 2001;41:347–366. doi: 10.1146/annurev.pharmtox.41.1.347. [PubMed] [Cross Ref]
  • Atkinson AJ, Colburn WA, DeGruttola VG, DeMets DL, Downing GJ, Hoth DF, Oates JA, Peck CC, Schooley RT, Spilker BA. et al. Biomarkers and surrogate endpoints: Preferred definitions and conceptual framework. Clinical Pharmacology & Therapeutics. 2001;69:89–95. [PubMed]
  • Gry M, Rimini R, Stromberg S, Asplund A, Ponten F, Uhlen M, Nilsson P. Correlations between RNA and protein expression profiles in 23 human cell lines. BMC Genomics. 2009;10:365. doi: 10.1186/1471-2164-10-365. [PMC free article] [PubMed] [Cross Ref]
  • Cho WC. Proteomics technologies and challenges. Genomics Proteomics Bioinformatics. 2007;5:77–85. doi: 10.1016/S1672-0229(07)60018-7. [PubMed] [Cross Ref]
  • Lubomirski M, D'Andrea MR, Belkowski SM, Cabrera J, Dixon JM, Amaratunga D. A consolidated approach to analyzing data from high-throughput protein microarrays with an application to immune response profiling in humans. J Comput Biol. 2007;14:350–359. doi: 10.1089/cmb.2006.0116. [PubMed] [Cross Ref]
  • Aebersold R, Mann M. Mass spectrometry-based proteomics. Nature. 2003;422:198–207. doi: 10.1038/nature01511. [PubMed] [Cross Ref]
  • Anderson NL. The roles of multiple proteomic platforms in a pipeline for new diagnostics. Mol Cell Proteomics. 2005;4:1441–1444. doi: 10.1074/mcp.I500001-MCP200. [PubMed] [Cross Ref]
  • Uhlen M, Oksvold P, Fagerberg L, Lundberg E, Jonasson K, Forsberg M, Zwahlen M, Kampf C, Wester K, Hober S. et al. Towards a knowledge-based Human Protein Atlas. Nat Biotechnol. 2010;28:1248–1250. doi: 10.1038/nbt1210-1248. [PubMed] [Cross Ref]
  • Kingsmore SF. Multiplexed protein measurement: technologies and applications of protein and antibody arrays. Nat Rev Drug Discov. 2006;5:310–320. doi: 10.1038/nrd2006. [PMC free article] [PubMed] [Cross Ref]
  • Schroder C, Jacob A, Tonack S, Radon TP, Sill M, Zucknick M, Ruffer S, Costello E, Neoptolemos JP, Crnogorac-Jurcevic T. et al. Dual-color proteomic profiling of complex samples with a microarray of 810 cancer-related antibodies. Mol Cell Proteomics. 2010;9:1271–1280. doi: 10.1074/mcp.M900419-MCP200. [PMC free article] [PubMed] [Cross Ref]
  • Borrebaeck CA, Wingren C. Recombinant antibodies for the generation of antibody arrays. Methods Mol Biol. 2011;785:247–262. doi: 10.1007/978-1-61779-286-1_17. [PubMed] [Cross Ref]
  • Gold L, Ayers D, Bertino J, Bock C, Bock A, Brody EN, Carter J, Dalby AB, Eaton BE, Fitzwater T. et al. Aptamer-based multiplexed proteomic technology for biomarker discovery. PLoS One. 2010;5:e15004. doi: 10.1371/journal.pone.0015004. [PMC free article] [PubMed] [Cross Ref]
  • Schwenk JM, Gry M, Rimini R, Uhlen M, Nilsson P. Antibody suspension bead arrays within serum proteomics. J Proteome Res. 2008;7:3168–3179. doi: 10.1021/pr700890b. [PubMed] [Cross Ref]
  • Boomsma D, Busjahn A, Peltonen L. Classical twin studies and beyond. Nat Rev Genet. 2002;3:872–882. doi: 10.1038/nrg932. [PubMed] [Cross Ref]
  • Vasan RS. Biomarkers of cardiovascular disease: molecular basis and practical considerations. Circulation. 2006;113:2335–2362. doi: 10.1161/CIRCULATIONAHA.104.482570. [PubMed] [Cross Ref]
  • Dieterle F, Ross A, Schlotterbeck G, Senn H. Probabilistic quotient normalization as robust method to account for dilution of complex biological mixtures. Application in 1H NMR metabonomics. Anal Chem. 2006;78:4281–4290. doi: 10.1021/ac051632c. [PubMed] [Cross Ref]
  • Wit E, McClure J. Statistics for Microarrays: Design, Analysis and Inference. John Wiley & Sons Australia, Limited; 2004.
  • Mahalanobis P. On the generalised distance in statistics. Proceedings of the National Institute of Sciences of India. 1936;2:49–55.
  • Stern AS, Magram J, Presky DH. Interleukin-12 an integral cytokine in the immune response. Life Sci. 1996;58:639–654. doi: 10.1016/S0024-3205(96)80003-8. [PubMed] [Cross Ref]
  • Porter S, Clark IM, Kevorkian L, Edwards DR. The ADAMTS metalloproteinases. Biochem J. 2005;386:15–27. doi: 10.1042/BJ20040424. [PubMed] [Cross Ref]
  • Goulter AB, Goddard MJ, Allen JC, Clark KL. ACE2 gene expression is up-regulated in the human failing heart. BMC Med. 2004;2:19. doi: 10.1186/1741-7015-2-19. [PMC free article] [PubMed] [Cross Ref]
  • Massague J, Andres J, Attisano L, Cheifetz S, Lopez-Casillas F, Ohtsuki M, Wrana JL. TGF-beta receptors. Mol Reprod Dev. 1992;32:99–104. doi: 10.1002/mrd.1080320204. [PubMed] [Cross Ref]
  • Portielje JEA, Kruit WHJ, Eerenberg AJM, Schuler M, Sparreboom A, Lamers CHJ, Bolhuis RLH, Stoter G, Huber C, Hack CE. Interleukin 12 induces activation of fibrinolysis and coagulation in humans. Brit J Haematol. 2001;112:499–505. doi: 10.1046/j.1365-2141.2001.02592.x. [PubMed] [Cross Ref]
  • Schwenk JM, Igel U, Kato BS, Nicholson G, Karpe F, Uhlen M, Nilsson P. Comparative protein profiling of serum and plasma using an antibody suspension bead array approach. Proteomics. 2010;10:532–540. doi: 10.1002/pmic.200900657. [PubMed] [Cross Ref]
  • Kricka LJ, Master SR. Quality control and protein microarrays. Clinical Chemistry. 2009;55:1053–1055. doi: 10.1373/clinchem.2009.126557. [PubMed] [Cross Ref]
  • Agaton C, Galli J, Hoiden Guthenberg I, Janzon L, Hansson M, Asplund A, Brundell E, Lindberg S, Ruthberg I, Wester K. et al. Affinity proteomics for systematic protein profiling of chromosome 21 gene products in human tissues. Mol Cell Proteomics. 2003;2:405–414. [PubMed]
  • Nilsson P, Paavilainen L, Larsson K, Odling J, Sundberg M, Andersson AC, Kampf C, Persson A, Al-Khalili Szigyarto C, Ottosson J. et al. Towards a human proteome atlas: high-throughput generation of mono-specific antibodies for tissue profiling. Proteomics. 2005;5:4327–4337. doi: 10.1002/pmic.200500072. [PubMed] [Cross Ref]
  • Colwill K, Persson H, Jarvik NE, Wyrzucki A, Wojcik J, Koide A, Kossiakoff AA, Koide S, Sidhu S, Dyson MR, A roadmap to generate renewable protein binders to the human proteome. Nat Methods. 2011.
  • Schwenk JM, Lindberg J, Sundberg M, Uhlen M, Nilsson P. Determination of Binding Specificities in Highly Multiplexed Bead-based Assays for Antibody Proteomics. Mol Cell Proteomics. 2007;6:125–132. [PubMed]
  • Hartley HO, Rao JKN. Maximum-likelihood estimation for the mixed analysis of variance model. Biometrika. 1967;54:93–108. [PubMed]
  • Searle SR, Casella G, McCulloch CE. Variance components. Wiley; 2006.
  • Rabe-Hesketh S, Skrondal A, Gjessing HK. Biometrical modeling of twin and family data using standard mixed model software. Biometrics. 2008;64:280–288. doi: 10.1111/j.1541-0420.2007.00803.x. [PubMed] [Cross Ref]
  • Bates D. Fitting linear mixed models in R. R News. 2005;5:27–30.
  • Nicholson G, Rantalainen M, Maher AD, Li JV, Malmodin D, Ahmadi KR, Faber JH, Hallgrimsdottir IB, Barrett A, Toft H. et al. Human metabolic profiles are stably controlled by genetic and environmental variation. Molecular systems biology. 2011;7:525. [PMC free article] [PubMed]
  • Akaike H. A new look at the statistical model identification. IEEE Transactions on Automatic Control. 1974;19:716–723. doi: 10.1109/TAC.1974.1100705. [Cross Ref]
  • Burnham KP, Anderson DR. Model selection and multimodel inference: a practical information-theoretic approach. 2. New York: Springer-Verlag; 2002.
Articles from Proteome Science are provided here courtesy of
BioMed Central