lgatto/synapter — pull request #68 by sgibb (2014-05-16 16:46:25+00:00)
Title: 3D grid search

Dear Laurent,

I am very sorry that I spent the whole day to implement the 3D grid search (#33). Please don't bite my head off!

###### Good news:
- The 3D grid search is implemented and the results of my test dataset are better than using the 2D grid search.
- Most methods are already adapted to the 3D grid search.

###### Bad news:
- `synergise` fails. `getGrid` returns an array instead of a matrix and I don't know how to put this array in the HTML report.
- `R CMD check` fails. `synapterTiny` was created with Synapter 0.9.0. I use `updateSynpterObject` in the adapted `searchGrid` and the grid search works like expected (`debug(synapterTiny$searchGrid` tells me that the grid results are correct) but the updated object is not returned and the `synapterTiny` object is not changed in the global environment. That's why all following functions fail, e.g. `getGrid` returns an empty list. A workaround would be to call `synapterTiny <- updateSynpaterObject(synapterTiny)` as first function. This problem occurs in all vignettes too. Do you know how we can call `updateSynapterObject` automatically and change the original object? (maybe with `<<-`, but I afraid of the sideeffects)

###### Changes:
- `searchGrid` has a new `imdiffs` argument (default: `seq(0.2, 2, 0.2)`)
- `plotFeatures` has a new `ionmobility` argument (default: `FALSE`) to plot three additional plots with *m/z* vs ion mobility.
- `plotGrid` has a new `maindim` argument to decide which dimension should be fixed in the heatmap respective for which dimension multiple heatmaps are created.
- `setImDiff`, `getImDiff` to set the difference threshold in ion mobility manually.
- `getBestGridParams`, `setBestGridParams`, `findEMRTs` are adapted to the new grid search.
- I also updated the documentation of these functions.

###### TODO:
- Fix all examples/vignettes which are based/using old synapter objects.
- Adapt `synergise`.
- `synergise` has not any `imdiff.from`, `imdiff.to` or `imdiff.by` argument, yet.
- Fix `synergise` HTML report for the array return value of `getGrid`.
- Find good default values for `imdiffs` in `searchGrid` and `setImDiff`.
- @pavel-shliaha has to test this feature for different datasets.
- Vignette?

Best wishes,

Sebastian

P.S.: The new grid search in action:
# gridSearch3

## Default Synapter Preprocessing

(same as in `synergise`)

```s
suppressPackageStartupMessages(library("synapter"))
packageVersion("synapter")
```

```
## [1] '1.7.4'
```

### File Import

```s
dataPath <- file.path("..", "..", "..", "playground", "data",
"biological_replicas_for_synapter")
iaPath <- file.path(dataPath, "HDMSE", "ia")
apexPath <- file.path(dataPath, "HDMSE", "apex")

input <- list(
identpeptide=file.path(iaPath,
"fermentor_02_sample_15_HDMSE_01_20140303114821",
"fermentor_02_sample_15_HDMSE_01_IA_final_peptide.csv"),
quantpeptide=file.path(iaPath,
"fermentor_03_sample_15_HDMSE_01_20140303120343",
"fermentor_03_sample_15_HDMSE_01_IA_final_peptide.csv"),
quantpep3d=file.path(apexPath,
"fermentor_03_sample_15_HDMSE_01_20140302183822",
"fermentor_03_sample_15_HDMSE_01_Pep3DAMRT.csv"),
fasta=file.path(dataPath,
"S.cerevisiae_uniprot_reference_canonical_16_09_13.fasta"))

## only for reproducibility
set.seed(1)
system.time(s <- Synapter(input))
```

```
## Reading identification final peptide file...
## Reading quantitation final peptide file...
## Reading quantitation Pep3D file...
## Filtering...
## Computing identification statistics...
```

```
## user system elapsed
## 115.311 1.396 117.252
```

### Filtering

```s
filterUniqueDbPeptides(s, verbose=TRUE)
```

```
## Cleavage results: 6621 proteins
## 179685 out of 298849 tryptic peptides are proteotypic.
## Peptide database file: ../../../playground/data/biological_replicas_for_synapter/S.cerevisiae_uniprot_reference_canonical_16_09_13.fasta
## Allowed mis-cleavages: 0
## I/L treatment: I != L
```

```s
filterPeptideLength(s, l = 7)

fdr <- fpr <- 0.01
ppm <- 20

filterIdentPepScore(s, method="BH", fdr=fdr)
filterQuantPepScore(s, method="BH", fdr=fdr)
filterIdentProtFpr(s, fpr=fpr)
filterQuantProtFpr(s, fpr=fpr)
filterIdentPpmError(s, ppm=ppm)
filterQuantPpmError(s, ppm=ppm)
```

### RT Model and Grid Search

```s
mergePeptides(s)

modelRt(s, span=0.05)

plotFeatures(s, what = "all", ionmobility = TRUE)
```
![rtmodel](https://cloud.githubusercontent.com/assets/1828443/2999680/fdad0c7c-dd18-11e3-8edd-4ba856e0a8f2.png)

```s

s
```

```
## Object of class "Synapter"
## Package version 1.7.4
## Data files:
## + Quantitation pep file: fermentor_03_sample_15_HDMSE_01_IA_final_peptide.csv
## + Identification pep file: fermentor_02_sample_15_HDMSE_01_IA_final_peptide.csv
## + Quantitation Pep3DAMRT file: fermentor_03_sample_15_HDMSE_01_Pep3DAMRT.csv
## + Fasta file: S.cerevisiae_uniprot_reference_canonical_16_09_13.fasta
## Log:
## [1] "Instance created on Fri May 16 16:58:16 2014"
## [2] "Read identification peptide data [32308,66]"
## [ 28 lines ]
## [31] "Merged identification and quantitation data [4150,140]"
## [32] "Modelled retention time using lowess and span 0.05"
```

```s
s2 <- s$copy()
s3 <- s$copy()
```

### Grid Search 2D (like the old one)

```s
system.time(searchGrid(s2, imdiffs = Inf, verbose=FALSE))
```

```
## user system elapsed
## 22.273 0.348 22.959
```

```s

getBestGridParams(s2)
```

```
## $prcntTotal
## nsd ppm imdiff
## [1,] 2 10 Inf
##
## $prcntModel
## nsd ppm imdiff
## [1,] 2 8 Inf
##
## $details
## nsd ppm imdiff
## [1,] 3 6 Inf
```

```s

setBestGridParams(s2)

findEMRTs(s2)

plotEMRTtable(s2)
```
![searchgrid2](https://cloud.githubusercontent.com/assets/1828443/2999685/0baf735a-dd19-11e3-8942-4901cda4365e.png)

```s
getEMRTtable(s2)
```

```
##
## 0 1 2 3 4 5 6 7 10
## 278 4235 519 65 8 1 3 1 1
```

```s
performance(s2)
```

```
## (S) Synapter: 4235 EMRTs uniquely matched.
## (I) Ident: 5111 peptides.
## (Q) Quant: 4974 peptides.
## Enrichment (S/Q): -14.86%
## Overlap:
## Q S QS
## 1492 753 3482
```

```s

s2
```

```
## Object of class "Synapter"
## Package version 1.7.4
## Data files:
## + Quantitation pep file: fermentor_03_sample_15_HDMSE_01_IA_final_peptide.csv
## + Identification pep file: fermentor_02_sample_15_HDMSE_01_IA_final_peptide.csv
## + Quantitation Pep3DAMRT file: fermentor_03_sample_15_HDMSE_01_Pep3DAMRT.csv
## + Fasta file: S.cerevisiae_uniprot_reference_canonical_16_09_13.fasta
## Log:
## [1] "Instance created on Fri May 16 16:58:16 2014"
## [2] "Read identification peptide data [32308,66]"
## [ 31 lines ]
## [34] "Set 'nsd', 'ppm error' and 'im diff' to 2, 8 and Inf (best 'auto')"
## [35] "Matched identification peptides and quantitation EMRTs [5111,127] (rescue)"
```

### Grid Search 3D

```s
system.time(searchGrid(s3, verbose=FALSE))
```

```
## user system elapsed
## 59.29 1.94 62.06
```

```s

getBestGridParams(s3)
```

```
## $prcntTotal
## nsd ppm imdiff
## [1,] 2 16 0.6
##
## $prcntModel
## nsd ppm imdiff
## [1,] 2 10 0.6
##
## $details
## nsd ppm imdiff
## [1,] 1 4 0.2
```

```s

setBestGridParams(s3)

plotGrid(s3, maindim = "im")
```
![searchgrid31](https://cloud.githubusercontent.com/assets/1828443/2999687/150c565c-dd19-11e3-91fa-a5cd082d9fc7.png)

```s
plotGrid(s3, maindim = "rt")
```
![searchgrid32](https://cloud.githubusercontent.com/assets/1828443/2999692/19f0a538-dd19-11e3-99b4-6a3ba8a4ffd2.png)

```s

findEMRTs(s3)

plotEMRTtable(s3)
```
![searchgrid33](https://cloud.githubusercontent.com/assets/1828443/2999694/1e630570-dd19-11e3-88e6-bd57cc43eba5.png)

```s
getEMRTtable(s3)
```

```
##
## 0 1 2 3 4
## 321 4678 102 9 1
```

```s
performance(s3)
```

```
## (S) Synapter: 4678 EMRTs uniquely matched.
## (I) Ident: 5111 peptides.
## (Q) Quant: 4974 peptides.
## Enrichment (S/Q): -5.95%
## Overlap:
## Q S QS
## 1132 836 3842
```

```s

s3
```

```
## Object of class "Synapter"
## Package version 1.7.4
## Data files:
## + Quantitation pep file: fermentor_03_sample_15_HDMSE_01_IA_final_peptide.csv
## + Identification pep file: fermentor_02_sample_15_HDMSE_01_IA_final_peptide.csv
## + Quantitation Pep3DAMRT file: fermentor_03_sample_15_HDMSE_01_Pep3DAMRT.csv
## + Fasta file: S.cerevisiae_uniprot_reference_canonical_16_09_13.fasta
## Log:
## [1] "Instance created on Fri May 16 16:58:16 2014"
## [2] "Read identification peptide data [32308,66]"
## [ 31 lines ]
## [34] "Set 'nsd', 'ppm error' and 'im diff' to 2, 10 and 0.6 (best 'auto')"
## [35] "Matched identification peptides and quantitation EMRTs [5111,127] (rescue)"
```