-
Notifications
You must be signed in to change notification settings - Fork 2
/
Copy pathms2dda.Rmd
181 lines (127 loc) · 6.68 KB
/
ms2dda.Rmd
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
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
---
title: "Processing DDA MS^2^ data - Basic Functions"
author: "Michael Witting"
date: "24th of June 2018"
output:
html_document: default
pdf_document: default
---
```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```
## Basic processing of DDA MS^2^ data with XCMS
Metabolite identification is usually achieved by comparing MS2 spectra against reference spectra. Such MS2 data can be generated using data dependent acquistion. For further processing it is important to link the MS^1^ features and collected MS^2^ spectra. Here simply means to isolate MS2 spectra belonging to a certain chromatographic peak are presented.
## Setup of XCMS
XCMS3 is able to use multiple cores (e.g. when working on multiple files). If required a parallel backend can be registered.
```{r, message=FALSE, warning=FALSE}
#load required library
library(xcms)
library(ggplot2)
#set number of cores
noOfCores <- 4
## Use socket based parallel processing on Windows systems
if (.Platform$OS.type == "unix") {
register(bpstart(MulticoreParam(noOfCores)))
} else {
register(bpstart(SnowParam(noOfCores)))
}
```
## Read MS data
The example file contains an RP-LC-MS run of Mix 1 taken from the Agilent Pesticide Reference mix measured on a Sciex TripleTof 6600. First we plot a base peak chromatogram to get an overview on the chromatogram.
```{r}
## MS1 and MS2 data and file has to be read only once.
xrms <- readMSData("data\\PestMix1_DDA.mzML", mode = "onDisk", centroided = TRUE)
#plot BPC
bpis <- chromatogram(xrms, aggregationFun = "max")
plot(bpis)
```
## Detect peaks
Chromatographic peak detection is performed using the CentWave algorithm. Since only a single sample is analyzed no retention time alignment is required. Peak detection is performed only on the MS^1^ level.
```{r}
#set parameters and find chromatographic peaks
ms1cwp <- CentWaveParam(snthresh = 5, noise = 100, ppm = 10, peakwidth = c(3,30))
ms1data <- findChromPeaks(xrms, param = ms1cwp, msLevel = 1)
#get all peaks
chromPeaks <- chromPeaks(ms1data)
#check detected peaks
plotChromPeaks(ms1data)
```
## Handling of DDA MS2 data
All MS^2^ spectra can be extracted by filtering out the MS^2^ level and retrieve all spectra using the <code>spectra()</code> function together with <code>filterMsLevel()</code>. A list of <code>Spectrum2</code> objects is returned (see [MSnbase Development](https://www.bioconductor.org/packages/devel/bioc/vignettes/MSnbase/inst/doc/MSnbase-development.html#27_spectrum_et_al:_classes_for_ms_spectra)).
```{r}
#get all MS2 spectra from DDA experiment
ms2spectra <- spectra(filterMsLevel(xrms, msLevel = 2))
#plot position of all acquired MS2 spectra
ms2pos <- data.frame(rt = unlist(lapply(ms2spectra, function(x) {return(x@rt)})),
mz = unlist(lapply(ms2spectra, function(x) {return(x@precursorMz)})))
plot(ms2pos$rt, ms2pos$mz)
```
MS^2^ spectra corresponding to a specific peak can be isolated by comparing the retention time of the peak with the acquisition timing of the MS^2^ spectrum. Additionally, the precursor m/z and peak m/z should match.
```{r consolidation}
#function for DDA data
#MS2 spectra need to be list of spectrum2 objects
getDdaMS2Scans <- function(chromPeak, ms2spectra, mzTol = 0.01, rtTol = 30) {
#make a new list
filteredMs2spectra <- list()
#isolate only the ones within range
for(i in 1:length(ms2spectra)) {
#isolate Ms2 spectrum
ms2spectrum <- ms2spectra[[i]]
#check if within range of peak
if(abs(chromPeak[4] - ms2spectrum@rt) < rtTol & abs(chromPeak[1] - ms2spectrum@precursorMz) < mzTol) {
filteredMs2spectra <- c(filteredMs2spectra, ms2spectrum)
}
}
#return list with filtered ms2 spectra
return(filteredMs2spectra)
}
```
Using this function MS^2^ spectra can be consolidated and only spectra that are within a certain limit around a detected peak are retained. This drastically reduces the number of spectra that have to be processed.
```{r}
#create empty list
filteredMs2spectra <- list()
# interate over all chromatographic peaks and get only spectra in range of peaks
for(i in 1:nrow(chromPeaks)) {
chromPeak <- chromPeaks[i,]
filteredMs2spectra_clipboard <- getDdaMS2Scans(chromPeak, ms2spectra)
filteredMs2spectra <- c(filteredMs2spectra, filteredMs2spectra_clipboard)
}
#plot the filtered spectra as red circles
ms2pos_filtered <- data.frame(rt = unlist(lapply(filteredMs2spectra, function(x) {return(x@rt)})),
mz = unlist(lapply(filteredMs2spectra, function(x) {return(x@precursorMz)})))
plot(ms2pos$rt, ms2pos$mz)
points(ms2pos_filtered$rt, ms2pos_filtered$mz, col = "red", cex = 2.5)
```
The <code>MSnbase</code> package provides functions for handling of MS^2^ spectra, e.g. creating a mirror plot of two spectra and to calculate their similarity.
```{r}
# make mirror plot of two spectra
plot(filteredMs2spectra[[100]], filteredMs2spectra[[101]])
# calculate dot product
compareSpectra(filteredMs2spectra[[100]], filteredMs2spectra[[101]], binSize = 0.01, fun = "dotproduct")
```
## Example Peak: Fluopicolide
Fluopicolide is used as example. An EIC for the [M+H]^+^ is generated using the chromatogram function. The function <code>getDdaMS2Scans</code> is used to isolate MS^2^ spectra that are belonging to Fluopicolide. For each MS^2^ collected for Fluopicolide a vertical red line is added to the plot to illustrate the timing of MS^2^ data collection.
```{r}
#isolate Fluopicolide
chromPeak <- chromPeaks[57,]
#Fluopicolide, exact mass = 381.965430576, [M+H]+ = 382.972706
eic <- chromatogram(xrms, aggregationFun = "max", mz = c(382.96, 382.98), rt = c(430,450))
plot(eic)
#filter out fitting spectra
filteredMs2spectra <- getDdaMS2Scans(chromPeak, ms2spectra)
#mark position in EIC
abline(v = unlist(lapply(filteredMs2spectra, function(x) {return(x@rt)})), col = "red")
```
Using functions from <code>MSnbase</code> similarity to a library spectrum can calculated. A example library spectrum for Fluopicolide was taken from the Metlin DB (MID 72270) ([Metlin](https://metlin.scripps.edu/)). In order to be used with the <code>compareSpectra</code> function it has to be stored as <code>Spectrum2</code> object. Afterwards it can be handled by <code>MSnbase</code>.
```{r}
#example spectrum (taken from Metlin MID 72270, ESI-Q-ToF)
librarySpectrum <- new("Spectrum2",
precursorMz = 382.9727,
mz = c(193.9949, 172.9555, 144.9603, 108.9841, 74.0161),
intensity = c(3, 100, 36, 10, 5),
centroided = TRUE)
# make mirror plot of two spectra
plot(filteredMs2spectra[[4]], librarySpectrum, xlim = c(50,400))
# calculate dot product
compareSpectra(filteredMs2spectra[[4]], librarySpectrum, binSize = 0.01, fun = "dotproduct")
```