Flagging False Positives Following Untargeted LC–MS

Nov 7, 2016 - Epigenetic changes can be studied with an untargeted characterization of histone post-translational modifications (PTMs) by liquid ...
0 downloads 0 Views 1MB Size
Subscriber access provided by University of Otago Library

Article

Flagging false positives following untargeted LCMS characterization of histone PTM combinations Sander Willems, Maarten Dhaenens, Elisabeth Govaert, Laura De Clerck, Paulien Meert, Christophe Van Neste, Filip Van Nieuwerburgh, and Dieter Deforce J. Proteome Res., Just Accepted Manuscript • DOI: 10.1021/acs.jproteome.6b00724 • Publication Date (Web): 07 Nov 2016 Downloaded from http://pubs.acs.org on November 12, 2016

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Journal of Proteome Research is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Proteome Research

Flagging untargeted

false

positives

LCMS

following

characterization

of

histone PTM combinations Sander Willems1,a; Maarten Dhaenens1,a; Elisabeth Govaerta; Laura De Clercka; Paulien Meerta; Christophe Van Nestea,b,c; Filip Van Nieuwerburgh2,a; Dieter Deforce2,a,* a:

Laboratory of Pharmaceutical Biotechnology, Ghent University, B-9000 Ghent, Belgium

b:

Bioinformatics Institute Ghent, Ghent University, 9052 Ghent, Belgium

c:

Center for Medical Genetics Ghent, Ghent University, 9000 Ghent, Belgium

1: Equal

contribution as first author

2:

Equal contribution as last author

*:

Corresponding author; Telephone: +32 (0)9 264 80 67; E-mail: [email protected]

1 ACS Paragon Plus Environment

Journal of Proteome Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Abstract Epigenetic changes can be studied with an untargeted characterization of histone posttranslational modifications (PTMs) by LCMS. While prior information about more than twenty types of histone PTMs exists, little is known about histone PTM combinations (PTMCs). Due to the combinatorial explosion it is intrinsically impossible to consider all potential PTMCs in a database search. Consequentially, high-scoring false positives with unconsidered but correct alternative isobaric PTMCs can occur. Current quality controls can neither estimate the amount of unconsidered alternatives nor flag potential false positives. Here, we propose a conceptual workflow that provides such options. In this workflow, an in silico modelling of all candidate isoforms with known-to-exist PTMs is made. The most frequently occurring PTM sets of these candidate isoforms are determined and used in several database searches. This suppresses the combinatorial explosion while considering as many candidate isoforms as possible. Finally, annotations can be classified as unique or ambiguous, the latter implying false positives. This workflow was evaluated on an LCMS dataset containing 44 histone extracts. We were able to consider 60% of all candidate isoforms. Importantly, 40% of all annotations were classified as ambiguous. This highlights the need for a more thorough evaluation of modified peptide annotations.

2 ACS Paragon Plus Environment

Page 2 of 40

Page 3 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Proteome Research

Keywords Combinatorial explosion; Histone; Mass spectrometry (MS); Post-translational modification (PTM); Quality control

3 ACS Paragon Plus Environment

Journal of Proteome Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Abbreviations and definitions cRAP

common Repository of Adventitious Proteins

D

Aspartic acid

Da

Dalton

E

Glutamic acid

FDR

False discovery rate

PTM

Post-translational modification

PTMC

Post-translational modification combination

K

Lysine

LC

Liquid chromatography

MS

Mass spectrometry

MSMS

Tandem mass spectrometry

m/z

Mass over charge ratio

PTM

Post-translational modification

R

Arginine

RT

Retention time

XIC

Extracted-ion chromatogram

4 ACS Paragon Plus Environment

Page 4 of 40

Page 5 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Proteome Research

Introduction Modifications of chromatin are essential in the regulation of many biological processes. A key element of this chromatin are the evolutionary conserved histone proteins1, which can get intensely modified with post-translational modifications (PTMs)2. Increasingly, the so-called histone code is becoming an important target in medical and pharmaceutical research3. The NIH Roadmap Epigenomics Mapping Consortium is therefore making efforts to systematically characterize the epigenomic landscape4. Within this landscape, some differential histone PTM expressions have already been used to indicate several types of cancer5. While still invaluable, immunoblotting assays are too targeted and thus inefficient for broader PTM characterizations of histone samples. For such characterizations, mass spectrometry (MS) coupled with liquid chromatography (LC) is gradually becoming the new standard6. Current MS acquisition workflows for histone analysis include histone extraction from a biological sample, chemical derivatization to block lysine residues, digestion of proteins into peptides using trypsin, a second round of propionylation to propionylate peptide N-termini, LC separation, and tandem MS (MSMS)7. Hereby, each run now yields data containing the LC retention times (RT), mass over charge ratios (m/z), and ion-fragment m/z of a few thousand different precursors. While this workflow is continuously being optimized8-13, awareness of potential pitfalls is essential for a thorough data analysis. Fully utilizing MS data remains a difficult task even in conventional proteomics studies, as evidenced by the high amount of unassigned MSMS spectra14. This amount is even higher for histones, as a result of the many different and frequent PTMs. These PTMs affect data acquisition by reduced spectrum quality due to e.g. reduced fragmentation of phosphorylated peptides15. As de novo and tag approaches primarily rely on consecutive b- and y-ions, and 5 ACS Paragon Plus Environment

Journal of Proteome Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

spectral histone libraries have yet to be made, database searches are still the prevalent approach16. Database searches make an in silico digestion of an anticipated protein database. They can allow for modifications by introducing variable or fixed amino acid modifications. Exhaustive lists of b- and y-ions are then created for each peptide, which in turn can be matched to experimental MSMS peaks in a spectrum for identification. However, each introduced variable modification increases the search-space exponentially for each single spectrum17. While error tolerant, eTag, and sequence tag approaches typically show which modifications are abundant in a sample, they often consider only a single modification per spectrum. Hereby these approaches lose information about PTM combinations (PTMCs), which have been hypothesized to be of significant influence in the histone code. Searching many spectra of hyper-modified histone peptides unfortunately implies that multiple variable modifications need to be searched simultaneously. This gives rise to a combinatorial explosion and consequentially unfeasibly long analysis times and computational demands. When characterizing a PTM-rich histone sample, a selection of PTMs has to be made. For an untargeted characterization however, there is no trivial selection of PTMs to consider, let alone a selection of PTMCs. Strategies that limit the number of investigated PTMs and accompanying combinatorial explosions have already been proposed18-20, but determining which PTMCs need to be searched mostly relies on the expertise of the analyst. While doing multiple searches and combining the results is more robust, this strategy in essence still suffers from educated guessing21. Even in conventional proteomics, some challenges remain after selection of search parameters. A well-known challenge is to determine the exact location of a single PTM. There 6 ACS Paragon Plus Environment

Page 6 of 40

Page 7 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Proteome Research

is often not enough information present in a spectrum to give a single correct answer. However, different solutions have been presented for this issue22, most notably D-, MD-, and PhosphoRS-scores23-25. As histones can be modified intensely, they present another related issue: the potential presence of many isobaric PTMCs. Many complex scenarios with isobaric PTMCs are possible. An example of such isobaric histone PTMCs are e.g. a combination of butyryl (C4H6O) and formyl (CO) with the same mass as a combination of propionyl (C3H4O) and acetyl (C2H2O). These modifications can be present on neighboring amino-acids and as a result the scores for both the butyryl+formyl and the propionyl+acetyl annotations differ due to only a few different b- or y-ions. Depending on search algorithms and noisy peaks, it is possible that the wrong annotation is assigned the higher score. Even worse, when only one of these two PTMCs is searched, the alternative correct annotation is not reported at all. This specific issue is rarely addressed and no solution has been presented to date. The amount of false positives in conventional proteomics is typically estimated by the false discovery rate (FDR). This FDR generally relies on target-decoy searches that only utilize peptide backbone sequences to estimate false positives. Therefore, FDRs are generally not adapted to handle hyper-modified peptides26. One approach with promising results that tries to tackle this issue is the global PTM identification of Shortreed, et al.27. Their main improvement over traditional database searching is the introduction of curated known-toexist PTMs. For histones in particular there are many known-to-exist PTMs2. Unfortunately Shortreed, et al.27 did not apply their approach on hyper-modified histone extracts, nor did they consider isobaric PTMCs. Furthermore this approach only calculates a global FDR and cannot flag individual false positives.

7 ACS Paragon Plus Environment

Journal of Proteome Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

In summary, histone analysis is computationally intensive and more susceptible to false positives than conventional proteomics. This is mainly due to the many different isobaric histone PTMCs. While solutions to some problems exist, there is no comprehensive histone data analysis workflow from start to end in which no bias is introduced by the used software and which estimates the quality of individual annotations. Here, we propose a conceptual workflow that provides such options for an untargeted characterization of histones. First, the amount of data is reduced by trimming spectra. Next, the size of the protein database is reduced by an initial database search. Hereafter, an in silico modelling at the MS level of all theoretical candidate isoforms with known-to-exist PTMs is made. This allows to determine the most frequently occurring PTMCs on these candidate isoforms, which are subsequently used in several database searches to identify experimental candidate isoforms at the MSMS level. Hereby the combinatorial explosion is reduced while considering as many candidate isoforms as possible. Finally, experimental candidate isoforms are assigned to precursors which are then individually classified to flag potential false positives. This conceptual workflow was used to evaluate a dataset containing 44 histone extracts with tools that best suit this particular dataset. However, the workflow was set up as general as possible and adjustments are easily plugged in. Examples of adjustments are for instance alternative sample preparation protocols with e.g. chemically induced methyl-esterification12 instead of amidation, other acquisition methodologies with e.g. electron-transfer dissociation instead of collision-induced dissociation, different software tools such as e.g. MaxQuant28 and MS Amanda29, and/or scripts other than Python such as e.g. C# and Java. Furthermore,

8 ACS Paragon Plus Environment

Page 8 of 40

Page 9 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Proteome Research

growing knowledge on known-to-exist PTMs is easy to introduce. For simplicity, these alternatives were not explored in this article.

9 ACS Paragon Plus Environment

Journal of Proteome Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 40

Material & Methods Dataset Dataset Exp44 from Govaert, et al.10 was used (ProteomeXchange ID: PXD002885). Briefly, data from 44 LCMS runs of a Jurkat cell line was obtained after: (a) extraction of histones by four different extraction protocols, (b) derivatization of the resulting protein sample with propionic anhydride to block free and mono-methylated lysines, (c) proteolytic digestion using trypsin, (d) a second round of derivatization to propionylate peptide N-termini, and (e) data acquisition of 400ng of each extract by high definition data dependent acquisition LCMS as described in Helm, et al.30 using a Synapt G2-Si. Effective propionylation of lysine (K) and peptide N-termini was validated as described in Meert, et al.8, thereby warranting Arg-C digestion specificity and fixed lysine and N-terminal propionylation. However, an undesired variable amidation on aspartic acid (D), glutamic acid (E), and the peptide C-termini was also found to be present due to the propionylation reaction8. Precursor mass tolerance was less than 10 ppm and fragment mass tolerance was 0.3 Dalton (Da) with predominantly b- and y-ions.

Data analysis workflow This dataset was analyzed with the following workflow: (i) feature detection and alignment, (ii) protein identification, (iii) candidate isoform modelling, (iv) PTM set prioritization, (v) candidate isoform identification, and (vi) feature classification (Figure 1). (i)

Feature detection and alignment

Raw data of all 44 runs was imported and aligned in Progenesis QIP 2.0 (Nonlinear Dynamics) with default settings. This resulted in several features, each described by a unique pair of 10 ACS Paragon Plus Environment

Page 11 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Proteome Research

precursor RT and precursor m/z, linked to a list of detected MSMS spectra. Features with less than three MSMS spectra were considered unusable for identification (see step (vi): Feature classification) and were ignored in subsequent analyses. Only features with charge between 2+ and 4+ and mass below 3,000 Da were retained to reduce subsequent computational resources. MSMS spectra were ranked according to their distance to the elution apex of the extracted-ion chromatogram (XIC), correlating to intensity and information content. The three spectra closest to the elution apex were retained for each feature and the remaining spectra were considered redundant. Peak-picking was done with the PLGS plugin in Progenesis QIP with default parameters. Finally, an *.mgf with all selected spectra as well as a list containing all retained feature masses were exported for further analysis. (ii)

Protein identification

The exported *.mgf from step (i) was searched with Mascot 2.5 (Matrix Science) and SearchGUI31 2.8.4 (MS-GF+32 and OMSSA33) to identify proteins which were co-extracted with histones. The database was set to SwissProt34 Human (downloaded April 25th, 2016) appended with the cRAP35database, containing common contaminating proteins in mass spectrometry, and a concatenated reversed decoy. Other parameters were set to: b- and y-ions, fixed propionylation of K and the peptide N-terminus, variable amidation of D, E, the peptide Cterminus and oxidation of methionine (M), precursor mass tolerance of 10 ppm, and fragment mass tolerance of 0.3 Da, Arg-C specificity with one allowed miscleavage, precursor charges between 2+ and 4+ and isotopes 0 and 1. The results of Mascot and SearchGUI were processed with PeptideShaker36 1.10.1 to obtain a list of confident proteins at a 1% FDR.

11 ACS Paragon Plus Environment

Journal of Proteome Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(iii)

Page 12 of 40

Candidate isoform modelling

A *.fasta file was created that contained all proteins from step (ii). This list was appended with all proteins in the cRAP35 database and all curated human histone variants present in SwissProt. Finally five Snapshot histones from the Cell Snapshot 2, a Cell paper entitled “SnapShot: Histone Modifications” containing many verified histone PTMs and their locations, were appended to form a complete protein database describing all present proteins. This *.fasta file was then updated to a *.peff (Proteomics Standards Initiative37 Extended Fasta Format) file where identities and locations of known-to-exist PTMs were included by appending them to sequence headers. Known-to-exist PTMs for Snapshot histones were parsed from the Cell Snapshot, while known-to-exist PTMs for all other proteins were obtained from SwissProt. However, histone variants have similar sequences that are redundant. This redundancy in the database was removed by aligning all SwissProt histones against histones from the Cell Snapshot2. Histones with more than 95% sequence identity against histones of the Cell Snapshot of Huang, et al.2 were collapsed into five unique Snapshot histone with additional point mutations defined as PTMs to retain sequence variants. Other histones with less similarity remained unchanged. Custom and/or newly discovered knownto-exist PTMs can easily be introduced in sequence headers in other experiments. An example of such a custom PTM is for instance a glycosylation with a unique fixed mass. Using an in-house Python package this *.peff file was in silico digested using Arg-C specificity with two allowed miscleavages. For each peptide backbone sequence, all candidate isoforms with all PTM combinations were generated. Variable chemical amidation of unmodified D, E and C-terminus as well as chemically fixed propionylation of unmodified and monomethylated K and N-terminus was added to each candidate isoform to form a search-space.

12 ACS Paragon Plus Environment

Page 13 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Proteome Research

Only candidate isoforms with mass below 3,000 Da and at most seven biological PTMs were retained to reduce subsequent computational resources. All selected feature masses from step (i) were mapped against this search-space, within a 10 ppm window, to obtain an exhaustive list of candidate isoforms at the MS level for each feature. (iv)

PTM set prioritization

The exported features from step (i) and their list of modelled candidate isoforms from step (iii) were used as input for the PTM set prioritization. Here, a search-space was created based on these candidate isoforms. While features were used to generate lists of candidate isoforms, they are not explicitly present within this search-space. It was assumed that all candidate proteins were equally likely for a feature, all peptide backbone sequences were equally likely for a protein, all PTM combinations were equally likely for a peptide backbone sequence, and all isoforms were equally likely for a PTM combination (Figure 2). Thus, each individual candidate isoform was given a weight equal to its relative frequency to form the final searchspace. More precisely, an in-house Python package was used to give each feature 𝐹𝑖 with at least a single candidate isoform weight 𝑤(𝐹𝑖 ), with the sum of all weights ∑ 𝑤(𝐹𝑖 ) equal to one. Within each feature, all candidate proteins 𝑃𝑖,𝑗 were given equal weight 𝑤(𝑃𝑖,𝑗 ) with the sum of these weights ∑ 𝑤(𝑃𝑖,𝑗 ) equal to the weight 𝑤(𝐹𝑖 ) of feature 𝐹𝑖 . This was iterated for each candidate peptide backbone sequence 𝑝𝑖,𝑗,𝑘 within each candidate protein 𝑃𝑖,𝑗 , each candidate PTM combination 𝑚𝑖,𝑗,𝑘,𝑙 within each candidate peptide backbone sequence 𝑝𝑖,𝑗,𝑘 , and each candidate isoform 𝐼𝑖,𝑗,𝑘,𝑙,𝑛 within each candidate PTM combination 𝑚𝑖,𝑗,𝑘,𝑙 .

13 ACS Paragon Plus Environment

Journal of Proteome Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

𝑤(𝐹𝑖 ) = 𝑤(𝑃𝑖,𝑗 ) =

𝑤(𝑝𝑖,𝑗,𝑘 ) =

Page 14 of 40

1 #𝐹𝑖 𝑤(𝐹𝑖 ) #(𝑃 in 𝐹𝑖 ) 𝑤(𝑃𝑖,𝑗 ) #(𝑝 in 𝑃𝑖,𝑗 )

𝑤(𝑚𝑖,𝑗,𝑘,𝑙 ) =

𝑤(𝑝𝑖,𝑗,𝑘 ) #(𝑚 in 𝑝𝑖,𝑗,𝑘 )

𝑤(𝐼𝑖,𝑗,𝑘,𝑙,𝑛 ) =

𝑤(𝑚𝑖,𝑗,𝑘.𝑙 ) #(𝐼 in 𝑚𝑖,𝑗,𝑘.𝑙 )

The result is a set of weights for all candidate isoforms 𝐼𝑖,𝑗,𝑘,𝑙,𝑛 that sums to one. A PTM set 𝑐 is now defined as a unique set of distinct types of PTMs. The size 𝑠(𝑐) of 𝑐 is equal to the amount of distinct PTMs. The weight 𝑤(𝑐) of 𝑐 equals the sum of weights ∑ 𝑤(𝑚𝑖,𝑗,𝑘,𝑙 ) of all candidate PTM combinations 𝑚𝑖,𝑗,𝑘,𝑙, with the same PTMs as 𝑐. The sum of weights of all PTM sets thereby equals one. ∑ 𝑤(𝑐𝑙 ) = ∑ ∑ 𝑤(𝑚𝑖,𝑗,𝑘,𝑙 ) = ∑ 𝑤(𝑚𝑖,𝑗,𝑘,𝑙 ) = 1 𝑙

𝑙 𝑖,𝑗,𝑘

𝑖,𝑗,𝑘,𝑙

When performing a database search for a PTM set 𝑐, all PTM subsets will equally be searched. The cumulative weight 𝑊 of a PTM set 𝑐 is therefore defined as the sum of all PTM sets that are a subset of 𝑐. 𝑊(𝑐) = ∑ 𝑤(𝑐𝑖 ) 𝑐𝑖 ⊂𝑐

Phosphorylations are known to generate poorly fragmented MSMS spectra following collisioninduced dissociation and are thus hard to identify15. Therefore, each cumulative weight 𝑊(𝑐) of a set c containing a phosphorylation was ignored in subsequent analyses. 14 ACS Paragon Plus Environment

Page 15 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Proteome Research

A greedy algorithm, thus locally optimal but not per se globally optimal, was then used to select multiple PTM sets with a maximum combined cumulative weight. To this end, PTM sets were ranked according to their cumulative weight. The set 𝑐best with highest cumulative weight was selected. To minimize overlap with other PTM sets, the weight 𝑤(𝑐𝑖 ) of all subsets 𝑐𝑖 contributing to the cumulative weight 𝑊(𝑐best ) were set to zero, and all cumulative weights of supersets containing 𝑐𝑖 were equally updated. Thereafter, the best PTM set with regards to this updated cumulative weight was selected and the cumulative weight of all its subsets was again set to zero before updating the cumulative weight of their supersets. This process of selecting optimal PTM sets and depleting subset weights was iterated 100 times. The resulting list with prioritized PTM sets was exported for subsequent sequential database searching. (v)

Candidate isoform identification

The exported *.mgf from step (i) was searched multiple times with Mascot and SearchGUI (MS-GF+ and OMSSA) with equal parameters, except the selected variable PTMs, which were chosen according to the exported list with prioritized PTM sets of step (iv). As database, the same *.fasta from step (iii) containing co-extracted proteins, all explicit human histone variants without collapsing homologs to unique Snapshot proteins, and cRAP proteins with a concatenated reverse decoy was used. Other parameters were set to: b- and y-ions, fixed propionylation of K and the peptide N-terminus, variable amidation of D, E and the peptide Cterminus, precursor mass tolerance of 10 ppm and fragment mass tolerance of 0.3 Da, Arg-C specificity with two allowed miscleavages, precursor charges between 2+ and 4+, and isotopes 0 and 1. The results were imported in different PeptideShaker projects and all candidate isoforms were validated at a 1% FDR. All confident candidate isoforms of all search engines were exported for all searches. As all these search engines relied on *.fasta files instead of

15 ACS Paragon Plus Environment

Journal of Proteome Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 40

*.peff files, they were unable to take prior PTM locations into account and some confident candidate isoforms can thus contain PTMs at previously unknown locations. (vi)

Feature classification

A classification was made for all features from step (i) by using the candidate isoforms from step (v). First, all candidate isoforms were assigned to their corresponding feature by m/z and RT. Here they were grouped at the protein level, peptide backbone sequence level and PTM combination level (Figure 2). For each isoform it was then determined whether it was reproducible, i.e. supported by at least two out of three spectra. Furthermore, all PTMs of each isoform were cross-validated with known-to-exist PTMs in SwissProt and the SnapShot. For each feature it was verified whether it had (A) a unique isoform, (B) a unique PTM combination or (C) a unique peptide sequence backbone. This verification was done using different subsets of isoforms: (1) only known and reproducible isoforms, (2) only known isoforms, (3) only reproducible isoforms, or (4) all isoforms. This resulted in 12 different feature classifications ranging from (1.A) having a unique known and reproducible isoform to (4.C) having a unique peptide backbone sequence. A final class (X.X) was used to indicate the remaining features that were never unique.

16 ACS Paragon Plus Environment

Page 17 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Proteome Research

Results & Discussion The histone dataset Exp44 from Govaert, et al.10 containing 44 LCMS runs from a Jurkat cell line was analyzed with (a) Progenesis QIP (Nonlinear Dynamics), (b) PeptideShaker 36, and (c) an in-house Python package according to the following workflow: (i) feature detection and alignment, (ii) protein identification, (iii) candidate isoform modelling, (iv) PTM set prioritization, (v) candidate isoform identification, and (vi) feature classification (Figure 1). Within this workflow steps (i) through (iii) primarily focus on data reduction and simplification. Step (iv) and (v) then maximize the feature coverage. Finally, the output is assessed qualitatively in step (vi).

Data reduction and simplification PTM-rich samples require much computational resources and long search times. As the data here originated from 44 extracts from the same sample, many spectra were redundant. Trimming these spectra thus allows to lower computational resources and search times. Furthermore, histone extracts consist of a limited protein composition and only a small protein database suffices for feature annotation. Finally, SwissProt and the Cell Snapshot2 contain much knowledge on potential modifications of these proteins. Modelling these modifications avoids considering rare or impossible feature annotations, again reducing computational resources and search times. (i)

Feature detection and alignment

Raw data of all runs was imported and aligned in Progenesis QIP to reduce data complexity. 166,708 MSMS spectra linked to 6,539 multiply charged features were detected and aligned over all samples. Features with larger masses and charges are increasingly difficult to annotate and model, so only features with mass below 3,000 Dalton and charge between 2+ and 4+ 17 ACS Paragon Plus Environment

Journal of Proteome Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 40

were retained for subsequent searches. Most of these features had multiple redundant spectra, so the MSMS spectra were trimmed. To ensure comparison between features, exactly three spectra were retained for each feature. The remaining 2,978 features and corresponding 8,934 (Figure 3) spectra were peak-picked and exported to an *.mgf for subsequent analysis (Supplementary File S-1). As the amount of spectra and search times are linearly correlated, a subsequent reduction in search-time of roughly 95% was achieved. For other experiments, alternatives such as e.g. MaxQuant can be used instead of Progenesis QIP to detect and align features and create an *.mgf. However, not all raw data can be processed with each tool. The raw data in this study was generated on a Waters Synapt G2-Si mass spectrometer and thus Progenesis QIP was used. This strategy also allowed optimal *.mgf peak-picking of Waters raw data with Waters’ PLGS plug-in in Progenesis QIP. (ii)

Protein identification

Proteins in the data runs were robustly identified with PeptideShaker by combining results from a variety of search algorithms. The *.mgf from step (i) was searched with Mascot and SearchGUI (OMSSA and MS-GF+) against all 20,247 proteins in the cRAP database and SwissProt Human. Only chemical modifications were used in the search parameters. Hereby 147 proteins were confidently identified at a 1% FDR in PeptideShaker (Supplementary Table S-1). This list was appended with the complete cRAP35 database and all known human histone variants present in SwissProt (Figure 3). All other proteins were considered not to be extracted from the original samples. This reduction of the protein database yields a reduction of 99% in terms of needed computational resources. Other search algorithm that accepts *.mgf as input could also be used, but this was not explored here.

18 ACS Paragon Plus Environment

Page 19 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Proteome Research

(iii)

Candidate isoform modelling

With this minimal protein database, an exhaustive list of theoretical candidate isoforms at the MS level was created. Twenty histones variants with more than 95% sequence identity with histones of the Cell Snapshot of Huang, et al.2 were collapsed into five Snapshot histone to remove redundancy. An in silico Arg-C digestion of the protein database was then made. This generated 11,083 peptide backbone sequences. 41 curated PTM types and their known-toexist locations were downloaded from SwissProt for proteins in the database. For the five snapshot proteins, 19 PTM types and their known-to-exist locations were parsed from the Cell Snapshot. Based on these PTM types and their known-to-exist locations combined with chemical modifications, 57,681,153 theoretical candidate isoforms with mass below 3,000 Dalton were exhaustively modelled (Supporting Information Figure S-1). If the known-to-exist locations of all PTM types were ignored, there was an estimated lower boundary of 109 for the number of theoretical candidate isoforms. For snapshot histone H4 this number of theoretical candidate isoforms was explicitly calculated to illustrate the difference between using known-to-exist locations and ignoring them. With known-to-exist locations, H4 was modelled to have 425,760 theoretical candidate isoforms using 19 different PTMs. Without these known-to-exist locations, as is the case when the modifications are set to variable on a certain amino acid, the number of theoretical candidate isoforms increased by more than a hundredfold to 65,728,215. Thus, while it is computationally challenging to calculate all theoretical candidate isoforms, using known-to-exist locations reduces the number of theoretical candidate isoforms by multiple orders of magnitude. It should be noted that the number of theoretical candidate isoforms with known-to-exist PTMs increases when novel PTMs are discovered and introduced.

19 ACS Paragon Plus Environment

Journal of Proteome Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 40

Finally, all feature masses were superimposed on all theoretical candidate isoforms. Hereby 7,055,399 different theoretical candidate isoforms were found to have the same mass within a 10 ppm window of a selected feature (Figure 3, Supporting Information Figure S-2, Supplementary File S-2). Of interest, there were 343 features without a single theoretical candidate isoform. Potential explanations are e.g. wrongfully filtered proteins, unknown PTMs, chemical noise, etcetera.

Maximizing feature coverage While steps (i) and (ii) effectively reduce the amount of spectra and minimize the protein database, step (iii) is only a theoretical simplification. Even so, after step (iii) there remain more candidate isoforms to search than feasible. It is therefore essential to formulate optimal search parameters that consider as many candidate isoforms as possible. This can be done by utilizing the theoretical candidate isoforms to prioritize PTM sets. Since a single search with optimal parameters still might not provide enough coverage of candidate isoforms, multiple sequential searches can be used. (iv)

PTM set prioritization

The relative frequency of all PTM sets was calculated to determine how many theoretical candidate isoforms each PTM set covers. The theoretical candidate isoforms of step (iii) had a total of 35,649 different PTM sets. However, not all PTM sets can efficiently be used in a database search, the main constraint being the size of a PTM set. SearchGUI advises not to use more than six variable PTMs to maintain manageable identification thresholds and search times. Propionylation of K and N-termini was considered as a fixed modifications and these modifications thus have no influence on the available PTM set size. Since there was chemical amidation of D, E and the peptide C-terminus, the maximum biological PTM set size was 20 ACS Paragon Plus Environment

Page 21 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Proteome Research

limited to three in this specific analysis. The total weight of all biological PTM sets of size four or larger equals 25% (Supplementary Table S-2). This 25% of theoretical candidate isoforms thus cannot be considered when only PTM sets of size three are considered, potentially yielding annotations with unconsidered but correct alternative isobaric PTMCs. When only biological PTM sets of size three or less are considered, there still remain 1,388 biological PTM sets. While multiple sequential searches can be used to consider more candidate isoforms, this number of biological PTM sets is still unfeasible to search. By selecting the top five of these biological PTM sets however, 60% of all theoretical candidate isoforms were already considered (Table 1). As the gain in candidate isoform coverage is limited when increasing this to six or more PTM sets, this was omitted here. Thus, the five following PTM sets were selected for subsequent analysis: (a) K-Acetyl, K-Methyl and K-Succinyl, (b) KTrimethyl, R-Methyl and R-Dimethyl, (c) K-Acetyl, K-Dimethyl and R-Citrulline, (d) K-Crotonyl, K-2-Hydroxyisobutyryl and R-Methyl, and (e) K-Formyl, K-Ubiquitin and S-Acetyl (Table 1). In every search, amidation on D, E and C-terminus were also set as variable to account for undesired chemical modifications. Phosphorylations were not prioritized in this study, as collision-induced dissociation was used for fragmentation and this makes phosphorylations hard to annotate. This constraint is not necessary for experiments containing spectra that are obtained with other fragmentation techniques. (v)

Candidate isoform identification

There is a plethora of search engines, each with a different scoring algorithm 38. Therefore, different search engines will have different outcomes, even though they have identical parameters. To compensate for such differences, SearchGUI (MS-GF+ and OMSSA) and Mascot were used in parallel to experimentally identify candidate isoforms in each spectrum

21 ACS Paragon Plus Environment

Journal of Proteome Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 40

at the MSMS level. The results were imported in PeptideShaker and validated at a 1% FDR level. The five PTM sets from step (iv) yielded 6,081 confident candidate isoforms (Figure 3, Supplementary Table S-3).

Quality assessment In step (v) the feature coverage was maximized in terms of considered candidate isoforms. However, quantity is not a measure of quality. There are for instance features that have several different candidate isoforms. This ambiguity can be present due to either a single spectrum annotated in a single search with a single search engine, or any combination of spectra, searches and engines. As search engines have different algorithms and different sequential searches were performed, there is no straightforward manner to compare scores of candidate isoforms. Even if this would have been possible, it is still ill-advised as MS spectra are noisy and contain much information that is not interpreted nor used. Furthermore, it is possible that some correct alternative candidate isoforms were never considered, leading to false positives scoring above significance threshold. It is therefore imperative to assess the quality of individual feature annotations by classifying them. A full overview of all features, their annotated isoforms, supporting spectra, supporting search engines, and supporting search files can be found in Supplementary Table S-3. (vi)

Feature classification

Candidate isoforms can be assigned to a feature at different levels: the peptide backbone sequence level, the PTM combination level, and the isoform level (Figure 2). The candidate isoforms identified in step (v) were assigned to 1,038 features (Supplementary Table S-4). In conventional proteomics the primary goal is often protein annotation, which requires confidently annotated peptides. For most annotated features this is not an issue, as evidenced 22 ACS Paragon Plus Environment

Page 23 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Proteome Research

by only 13% of features with ambiguous peptide backbone sequence annotations (Figure 4A). It should be noted that a small protein database was used and thus the search-space contained only a few peptides. Other common proteomics studies try to determine the exact location of PTMs that are e.g. regulators in several pathways. This has proven to be difficult in our samples, since only 21% of all features have exactly one annotated isoform (Figure 4A). Several tools have already been developed to address the specific issue of determining exact PTM locations and are discussed elsewhere23-25. Here however, the exact location of the PTMs is not the main cause for ambiguous isoforms. The 5% unique PTMC annotations (Figure 4A) that do not have a unique isoform, imply that localization is not the major issue. A major cause for ambiguity at the isoform level seems to be at the PTM combination level, as more than 60% of all features are annotated with a unique peptide backbone sequence, but with different PTMCs (Figure 4A). These results indicate that there are many isobaric PTMCs, within 10 ppm. This rarely-discussed issue is presumably more abundant in histone samples as these are subject to many PTMs. While it is possible that features contain chimeric spectra leading to ambiguity, it is more likely that some candidate isoforms are false positives. Chimeric features can potentially be distinguished from false positives by an increased RT width (Supporting Information Figure S-4), but this falls outside the scope of this article. Pre-filtering irreproducible candidate isoforms To tackle this ambiguity, two criteria were considered to pre-filter false positive isoforms before feature classification. The first criterion assumes that false positives are more prone to chance and thus less reproducible than true positives. All candidate isoforms that are supported by only one spectrum were thus considered as potential false positives and removed. This criterion is similar to conventional protein inference, where two or more peptides are required to infer presence of a protein. This criterion is justifiable in this particular 23 ACS Paragon Plus Environment

Journal of Proteome Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 40

study as 44 runs were merged and reproducibility was expected for the majority of features. Note that with exactly three spectra per feature and minimum two supporting spectra per isoform, feature ambiguity implies that at least one spectrum is also ambiguous and can be annotated with multiple isoforms. When classifying features after filtering irreproducible candidate isoforms (Supporting Information Figure S-3), relatively less features were found to be ambiguous at the peptide backbone sequence level (Figure 4B). However, ambiguity at the PTM combination and the isoform level remained the same (Figure 4B). Alternative definitions of reproducibility, e.g. requiring three supporting spectra instead of two, have similar results (Supporting Information Figure S-5). Pre-filtering unknown candidate isoforms Known PTM locations were used as a second independent criterion to filter false positive isoforms. While candidate isoforms with known-to-exist PTM locations were modelled in step (iii), the used search engines lose such information when PTMs are defined as variable. Therefore, some candidate isoforms were experimentally identified with PTMs at unknown locations. While it is possible that novel PTM locations were indeed identified, this requires a thorough validation39. This second criterion is similar to the use of decoy databases, where reversed or randomized peptide sequence are assumed to be false positives as they are a priori unknown. Removing candidate isoforms that are unknown solves almost all ambiguity at the peptide backbone sequence level (Figure 4C). Furthermore, the percentage of features with unique isoforms almost doubles by filtering unknown isoforms and the amount of features with a unique PTMCs more than triples (Figure 4C). This results in more than 55% of all annotated features to be unique at either the PTM combination level or the isoform level, as opposed to less than 25% when no filtering is applied. Equally important is that 927 features out of the original 1038 features still retain at least a single isoform after filtering out unknown 24 ACS Paragon Plus Environment

Page 25 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Proteome Research

isoforms (Figure 4C). Thus the improvements are not only relative but also absolute. It should be noted that the third search of step (v) includes chemical amidation of D, E and the Cterminus as well as biological citrullination of R. As these PTMs have opposite masses of ±0.9840, they have a net mass of zero when both present. Therefore, search engines can often annotate unmodified peptides with similar score as their counterpart with both amidation and citrullination. Since citrullination is a biological PTM, isoforms with this PTM can be filtered out if they are unknown (Supplementary Table S-3). Including this filtering step in our workflow thus makes the data analysis more robust to suboptimal sample preparations. Pre-filtering irreproducible and unknown candidate isoforms In a fourth and final classification, only isoforms that were both reproducible and known were retained. This resulted in only 1% of all features to be ambiguous at the peptide backbone sequence level (Figure 4D, Supporting Information Figure S-5), while almost 60% of all annotated features were classified as being unique at the PTM combination level or even isoform level (Figure 4D, Supporting Information Figure S-5). With respect to the unfiltered results, this corresponds to more than a twofold reduction in ambiguity.

25 ACS Paragon Plus Environment

Journal of Proteome Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 40

Conclusion Histone post-translational modifications are a key element in epigenetics. However, an untargeted characterization of histone extracts using a bottom-up MS approach is computationally intensive and susceptible to false positives. With our conceptual workflow we were able to reduce the needed computational resources several orders of magnitude by trimming spectra, filtering proteins, and modelling candidate isoforms. However, even with this reduction it is intrinsically impossible to consider all candidate isoforms. In contrast to other approaches, our workflow enables the user to minimize, estimate and report the amount of unconsidered candidate isoforms. In a final classification step, 40% of all annotated features were ambiguously annotated with isobaric PTMCs and thus flagged as potential false positives. We conclude that there are many features that can be annotated with isobaric PTMCs. Our workflow is successfully able to detect such features and flag their annotations as potential false positives. Moreover, the high amount of potential false positives indicates that such a workflow is imperative for a correct interpretation of histone PTM characterizations.

26 ACS Paragon Plus Environment

Page 27 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Proteome Research

References (1)

Jenuwein, T.; Allis, C. D. Translating the Histone Code. Science 2001, 293, (5532), 1074-1080.

(2)

Huang, H.; Sabari, B. R.; Garcia, B. A.; Allis, C. D.; Zhao, Y. SnapShot: Histone Modifications. Cell 2014, 159, (2), 458-458.e1.

(3)

Priestley, C. C.; Anderton, M.; Doherty, A. T.; Duffy, P.; Mellor, H. R.; Powell, H.; Roberts, R. Epigenetics - relevance to drug safety science. Toxicology Research 2012, 1, (1), 23-31.

(4)

Bernstein, B. E.; Stamatoyannopoulos, J. A.; Costello, J. F.; Ren, B.; Milosavljevic, A.; Meissner, A.; Kellis, M.; Marra, M. A.; Beaudet, A. L.; Ecker, J. R.; Farnham, P. J.; Hirst, M.; Lander, E. S.; Mikkelsen, T. S.; Thomson, J. A. The NIH Roadmap Epigenomics Mapping Consortium. Nature biotechnology 2010, 28, (10), 1045-1048.

(5)

LeRoy, G.; DiMaggio, P.; Chan, E.; Zee, B.; Blanco, M.; Bryant, B.; Flaniken, I.; Liu, S.; Kang, Y.; Trojer, P.; Garcia, B. A quantitative atlas of histone modification signatures from human cancer cells. Epigenetics & Chromatin 2013, 6, (1), 20.

(6)

Britton, L.-M. P.; Gonzales-Cope, M.; Zee, B. M.; Garcia, B. A. Breaking the histone code with quantitative mass spectrometry. Expert review of proteomics 2011, 8, (5), 631-643.

(7)

Huang, H.; Lin, S.; Garcia, B. A.; Zhao, Y. Quantitative Proteomic Analysis of Histone Modifications. Chemical Reviews 2015, 115, (6), 2376-2418.

(8)

Meert, P.; Govaert, E.; Scheerlinck, E.; Dhaenens, M.; Deforce, D. Pitfalls in histone propionylation during bottom-up mass spectrometry analysis. PROTEOMICS 2015, 15, (17), 2966–2971.

27 ACS Paragon Plus Environment

Journal of Proteome Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(9)

Contrepois, K.; Ezan, E.; Mann, C.; Fenaille, F. Ultra-High Performance Liquid Chromatography−Mass Spectrometry for the Fast Profiling of Histone PostTranslational Modifications. Journal of Proteome Research 2010, 9, (10), 5501-5509.

(10)

Govaert, E.; Van Steendam, K.; Scheerlinck, E.; Vossaert, L.; Meert, P.; Stella, M.; Willems, S.; De Clerck, L.; Dhaenens, M.; Deforce, D. Extracting histones for the specific purpose of label-free MS. PROTEOMICS in press, 10.1002/pmic.201600341.

(11)

Meert, P.; Dierickx, S.; Govaert, E.; De Clerck, L.; Willems, S.; Dhaenens, M.; Deforce, D. Tackling aspecific side reactions during histone propionylation: The promise of reversing overpropionylation. PROTEOMICS 2016, 16, (14), 1970-1974.

(12)

Paternoster, V.; Edhager, A. V.; Sibbersen, C.; Nielsen, A. L.; Børglum, A. D.; Christensen, J. H.; Palmfeldt, J. Quantitative assessment of methyl-esterification and other side reactions in a standard propionylation protocol for detection of histone modifications. PROTEOMICS 2016, 16, (14), 2059-2063.

(13)

Soldi, M.; Cuomo, A.; Bonaldi, T. Quantitative assessment of chemical artefacts produced by propionylation of histones prior to mass spectrometry analysis. PROTEOMICS 2016, 16, (14), 1952-1954.

(14)

Ning, K.; Fermin, D.; Nesvizhskii, A. I. Computational analysis of unassigned highquality MS/MS spectra in proteomic data sets. PROTEOMICS 2010, 10, (14), 27122718.

(15)

Boersema, P. J.; Mohammed, S.; Heck, A. J. R. Phosphopeptide fragmentation and analysis by mass spectrometry. Journal of Mass Spectrometry 2009, 44, (6), 861-878.

(16)

Sidoli, S.; Cheng, L.; Jensen, O. N. Proteomics in chromatin biology and epigenetics: Elucidation of post-translational modifications of histone proteins by mass spectrometry. Journal of Proteomics 2012, 75, (12), 3419-3433. 28 ACS Paragon Plus Environment

Page 28 of 40

Page 29 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(17)

Journal of Proteome Research

Schwämmle, V.; Verano-Braga, T.; Roepstorff, P. Computational and statistical methods for high-throughput analysis of post-translational modifications of proteins. Journal of Proteomics 2015, 129, 3-15.

(18)

Bern, M.; Cai, Y.; Goldberg, D. Lookup Peaks:  A Hybrid of de Novo Sequencing and Database Search for Protein Identification by Tandem Mass Spectrometry. Analytical Chemistry 2007, 79, (4), 1393-1400.

(19)

Na, S.; Paek, E. Software eyes for protein post-translational modifications. Mass Spectrometry Reviews 2015, 34, (2), 133-147.

(20)

Kertész-Farkas, A.; Reiz, B.; Vera, R.; Myers, M. P.; Pongor, S. PTMTreeSearch: a novel two-stage tree-search algorithm with pruning rules for the identification of posttranslational modification of proteins in MS/MS spectra. Bioinformatics 2014, 30, (2), 234-241.

(21)

Yuan, Z.-F.; Lin, S.; Molden, R. C.; Garcia, B. A. Evaluation of Proteomic Search Engines for the Analysis of Histone Modifications. Journal of Proteome Research 2014, 13, (10), 4470-4478.

(22)

Chalkley, R. J.; Clauser, K. R. Modification Site Localization Scoring: Strategies and Performance. Molecular & Cellular Proteomics 2012, 11, (5), 3-14.

(23)

Vaudel, M.; Breiter, D.; Beck, F.; Rahnenführer, J.; Martens, L.; Zahedi, R. P. D-score: A search engine independent MD-score. PROTEOMICS 2013, 13, (6), 1036-1041.

(24)

Savitski, M. M.; Lemeer, S.; Boesche, M.; Lang, M.; Mathieson, T.; Bantscheff, M.; Kuster, B. Confident Phosphorylation Site Localization Using the Mascot Delta Score. Molecular & Cellular Proteomics 2011, 10, (2), 10.1074/mcp.M110.003830.

29 ACS Paragon Plus Environment

Journal of Proteome Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(25)

Taus, T.; Köcher, T.; Pichler, P.; Paschke, C.; Schmidt, A.; Henrich, C.; Mechtler, K. Universal and Confident Phosphorylation Site Localization Using phosphoRS. Journal of Proteome Research 2011, 10, (12), 5354-5362.

(26)

Jeong, K.; Kim, S.; Bandeira, N. False discovery rates in spectral identification. BMC Bioinformatics 2012, 13, (16), 1-15.

(27)

Shortreed, M. R.; Wenger, C. D.; Frey, B. L.; Sheynkman, G. M.; Scalf, M.; Keller, M. P.; Attie, A. D.; Smith, L. M. Global Identification of Protein Post-translational Modifications in a Single-Pass Database Search. Journal of Proteome Research 2015, 14, (11), 4714-4720.

(28)

Cox, J.; Mann, M. MaxQuant enables high peptide identification rates, individualized p.p.b.-range mass accuracies and proteome-wide protein quantification. Nat Biotech 2008, 26, (12), 1367-1372.

(29)

Dorfer, V.; Pichler, P.; Stranzl, T.; Stadlmann, J.; Taus, T.; Winkler, S.; Mechtler, K. MS Amanda, a Universal Identification Algorithm Optimized for High Accuracy Tandem Mass Spectra. Journal of Proteome Research 2014, 13, (8), 3679-3684.

(30)

Helm, D.; Vissers, J. P. C.; Hughes, C. J.; Hahne, H.; Ruprecht, B.; Pachl, F.; Grzyb, A.; Richardson, K.; Wildgoose, J.; Maier, S. K.; Marx, H.; Wilhelm, M.; Becher, I.; Lemeer, S.; Bantscheff, M.; Langridge, J. I.; Kuster, B. Ion Mobility Tandem Mass Spectrometry Enhances Performance of Bottom-up Proteomics. Molecular & Cellular Proteomics 2014, 13, (12), 3709-3715.

(31)

Vaudel, M.; Barsnes, H.; Berven, F. S.; Sickmann, A.; Martens, L. SearchGUI: An opensource graphical user interface for simultaneous OMSSA and X!Tandem searches. PROTEOMICS 2011, 11, (5), 996-999.

30 ACS Paragon Plus Environment

Page 30 of 40

Page 31 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(32)

Journal of Proteome Research

Kim, S.; Pevzner, P. A. MS-GF+ makes progress towards a universal database search tool for proteomics. Nature Communications 2014, 5, 10.1038/ncomms6277.

(33)

Geer, L. Y.; Markey, S. P.; Kowalak, J. A.; Wagner, L.; Xu, M.; Maynard, D. M.; Yang, X.; Shi, W.; Bryant, S. H. Open Mass Spectrometry Search Algorithm. Journal of Proteome Research 2004, 3, (5), 958-964.

(34)

The UniProt, C. UniProt: a hub for protein information. Nucleic Acids Research 2015, 43, (D1), D204-D212.

(35)

Craig, R.; Cortens, J. P.; Beavis, R. C. Open Source System for Analyzing, Validating, and Storing Protein Identification Data. Journal of Proteome Research 2004, 3, (6), 1234-1242.

(36)

Vaudel, M.; Burkhart, J. M.; Zahedi, R. P.; Oveland, E.; Berven, F. S.; Sickmann, A.; Martens, L.; Barsnes, H. PeptideShaker enables reanalysis of MS-derived proteomics data sets. Nat Biotech 2015, 33, (1), 22-24.

(37)

Orchard, S.; Hermjakob, H.; Apweiler, R. The Proteomics Standards Initiative. PROTEOMICS 2003, 3, (7), 1374-1376.

(38)

Nesvizhskii, A. I. A survey of computational methods and error rate estimation procedures for peptide and protein identification in shotgun proteomics. Journal of Proteomics 2010, 73, (11), 2092-2123.

(39)

Lee, S.; Tan, M.; Dai, L.; Kwon, O. K.; Yang, J. S.; Zhao, Y.; Chen, Y. MS/MS of Synthetic Peptide Is Not Sufficient to Confirm New Types of Protein Modifications. Journal of Proteome Research 2013, 12, (2), 1007-1013.

31 ACS Paragon Plus Environment

Journal of Proteome Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 32 of 40

Acknowledgements This research was primarily funded by the Research Foundation Flanders (FWO) via research project grant G013916N. Partial funding was received through FWO mandate 12E9716N awarded to Maarten Dhaenens, Innovation through Science and Technology in Flanders (IWT) mandate SB-11179 awarded to Paulien Meert, IWT mandate SB-141209 awarded to Laura De Clerck, and Special Research Fund (BOF) mandate 01D23313 awarded to Elisabeth Govaert. The authors would like to express their gratitude to Marc Vaudel for his support on PeptideShaker issues and to Lennart Martens for proofreading the manuscript.

32 ACS Paragon Plus Environment

Page 33 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Proteome Research

Conflict of interest The authors declare no competing financial interest.

33 ACS Paragon Plus Environment

Journal of Proteome Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 34 of 40

Supporting information The following files are available free of charge at ACS website http://pubs.acs.org:        

Supporting_Information.pdf Supplementary_Table_S-1_Protein_Identification.txt Supplementary_Table_S-2_PTM_Set_Prioritization.txt Supplementary_Table_S-3_Experimental_Candidate_Isoforms.txt Supplementary_Table_S-4_Classified_Features.txt Supplementary_File_S-1_Progenesis_Spectra.mgf.txt.zip Supplementary_File-S-2_Theoretical_Candidate_Isoforms.xml.txt.zip Supplementary_File-S-3_Source_Code.zip

Additional information on supplementary tables and files can be found in Supporting Information.

34 ACS Paragon Plus Environment

Page 35 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Proteome Research

Tables and figures Figures Figure 1: Workflow

An LCMS dataset containing 44 LCMS runs of different histone extracts of a Jurkat cell line was analyzed with the following workflow: (i) Raw data is first simplified by aligning all precursors and reduced by trimming MSMS spectra to exactly three per feature. (ii) These spectra are then used to obtain a minimal database representing all proteins in the extract. (iii) This database is subsequently combined with known-to-exist PTMs to exhaustively model all candidate isoforms at the MS level. (iv) PTM combinations on these candidate isoforms are then weighted to prioritize PTM sets for subsequent database searches, thereby suppressing a combinatorial explosion. (v) Multiple PTM sets are sequentially searched to annotate spectra while considering as many candidate isoforms as possible. (vi) Features are classified to determine whether they are uniquely annotated at the peptide level, PTM combination level or isoform level when requiring candidate isoform to be reproducible or known. For this 35 ACS Paragon Plus Environment

Journal of Proteome Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 36 of 40

specific dataset the following software tools were used: (a) the commercial MS analysis software tool Progenesis QIP 2.0 (Nonlinear Dynamics) as alignment software, (b) the opensource Java package PeptideShaker36 in combination with SearchGUI31 and Mascot (Matrix Science) as MS2 identification software, and (c) Python as custom programming language. Figure 2: Feature organization

Features can be organized in a tree-based manner, depending on their candidate isoforms. Each feature can thereby encompass different proteins. For each protein different peptide backbone sequences can be annotated. Each peptide backbone sequence can be modified with different PTM combinations. Finally, different isoforms describe different locations of a PTM combination on a peptide. A PTM combination is always coupled to a specific peptide backbone sequence, while a PTM set is independent of a peptide.

36 ACS Paragon Plus Environment

Page 37 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Proteome Research

Figure 3: Filtering effects

The histone analysis workflow includes six steps, of which the first three are primarily focused on reducing and simplifying data. Hereby much redundant or unnecessary data is discarded. The data output of each step is used as data input in subsequent steps of the workflow. Note that the discarded number of theoretical candidate isoforms is not explicitly calculated as this would be too computationally intensive. Instead, a conservative estimate was made.

37 ACS Paragon Plus Environment

Journal of Proteome Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 38 of 40

Figure 4: Feature classification

All features are classified to their level of uniqueness. A feature is unique at the isoform level if it has only one candidate isoform over three MSMS spectra. A feature is unique at the PTM combination level if all its candidate isoforms differ only in the location of their PTMs. If a feature has candidate isoforms with different PTM combinations but all on the same peptide backbone sequence, it is unique at the peptide backbone sequence level. All other features are classified as not being unique. Each candidate isoform itself can be categorized as reproducible when it is supported by at least two different spectra. A candidate isoform can also be categorized as known, when the locations of all its PTMs are known in SwissProt 34 or the Snapshot of Huang, et al.2. The feature classification was made for all features with at least a single candidate isoform when (A) all isoforms are retained, (B) only reproducible isoforms

38 ACS Paragon Plus Environment

Page 39 of 40

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Proteome Research

were retained, (C) only known isoforms were retained, and (D) only isoforms that were both reproducible and known were retained. For TOC only

Tables Table 1: Prioritized PTM combinations Priority 1 2 3 4 5 6 7 8 9 10

PTM combination K-Ac K-Me3 K-Ac K-Cr K-Fo K-Ac K-Ac K-Su K-Me K-Ac

K-Me R-Me K-Me2 K-Hib K-Ub K-Me K-Cr R-Me K-Me3 R-Me2

K-Su R-Me2 R-Cit R-Me S-Ac R-Me R-Me2 Y-OH R-Cit S-Ac

Cumulative Weight (%) Original Depleted Sequential 55.51 55.51 55.51 51.87 2.06 57.58 54.45 1.22 58.80 51.37 1.02 59.82 50.59 0.78 60.60 55.42 0.49 61.09 54.51 0.48 61.57 51.59 0.36 61.94 51.42 0.35 62.28 54.53 0.35 62.63

39 ACS Paragon Plus Environment

Journal of Proteome Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

11 12 13 14 15 16 17 18 19 20

K-Me K-Cr K-Bu K-Hib K-Cr K-Fo K-Ac K-Hib K-Cr K-Ac

K-Me2 K-Fo K-Ub K-Me2 K-Me K-Pr K-Hib K-Ub K-Me3 K-Cr

50

K-Bu

K-Me3

100

K-Ac

K-Me

R-Me2 K-Me2 R-Me R-Me T-Ac R-Me K-Me3 R-Me2 R-Me2 R-Cit … R-Me … S-Ac

Page 40 of 40

51.86 50.99 50.86 51.47 51.13 50.81 53.89 50.94 51.27 54.36

0.33 0.32 0.31 0.27 0.25 0.23 0.21 0.19 0.19 0.18

62.96 63.27 63.58 63.85 64.10 64.33 64.54 64.73 64.92 65.10

51.18

0.05

67.57

54.84

0.03

69.33

All candidate isoforms with known-to-exist post-translational modifications (PTMs) were exhaustively modelled for all features. All PTM sets were then weighted based on these candidate isoforms (Materials & Methods). The cumulative weight expresses the summed weight of a PTM set and all its subsets, including the empty set (Original). Subsequently, a greedy algorithm is used to prioritize PTM sets so that their coverage is maximal when combined with each other (Sequential). This algorithm selects the PTM combination with highest weight and the weight of all its subsets is depleted to zero. The weight of each PTM set is then recalculated (Depleted). A second PTM combination with highest depleted weight is selected and all PTM subsets are again depleted. This process is iterated several times. Note that the first PTM combination is the only PTM combination where the empty set has not yet been depleted, making it excessively large compared to subsequent PTM sets.

40 ACS Paragon Plus Environment