The following content is descibed in more detail in @NestLink (under review NMETH-A35040).

library(NestLink)
stopifnot(require(specL))

Input Pool Frequency

aa_pool_x8 <- c(rep('A', 12), rep('S', 0), rep('T', 12), rep('N', 12),
    rep('Q', 12), rep('D', 8),  rep('E', 0), rep('V', 12), rep('L', 0),
    rep('F', 0), rep('Y', 8), rep('W', 0), rep('G', 12), rep('P', 12))

aa_pool_1_2_9_10 <- c(rep('A', 8), rep('S', 7), rep('T', 7), rep('N', 6),
    rep('Q', 6), rep('D', 8), rep('E', 8), rep('V', 9), rep('L', 6),
    rep('F', 5), rep('Y', 9), rep('W', 6),  rep('G', 15), rep('P', 0))

aa_pool_3_8 <- c(rep('A', 5), rep('S', 4), rep('T', 5), rep('N', 2),
    rep('Q', 2), rep('D', 8), rep('E', 8), rep('V', 7), rep('L', 5),
    rep('F', 4), rep('Y', 6), rep('W', 4),  rep('G', 12), rep('P', 28))

Sanity Check

table(aa_pool_x8)
length(aa_pool_x8)
table(aa_pool_1_2_9_10)
length(aa_pool_1_2_9_10)
table(aa_pool_3_8)
length(aa_pool_3_8)

Compose Peptides

GPGXXXXXXXX(VR|VSR|VFGIR|VSGER) peptide

replicate(10, compose_GPGx8cTerm(pool=aa_pool_x8))

GPYYXXXXXXYYR peptide

compose_GPx10R(aa_pool_1_2_9_10, aa_pool_3_8)

Generate peptides

set.seed(2)
(sample.size <- 3E+04)
peptides.GPGx8cTerm <- replicate(sample.size, compose_GPGx8cTerm(pool=aa_pool_x8))
peptides.GPx10R <- replicate(sample.size, compose_GPx10R(aa_pool_1_2_9_10, aa_pool_3_8))
# write.table(peptides.GPGx8cTerm, file='/tmp/pp.txt')

Peptide mass

Compute peptide mass

library(protViz)
(smp.peptide <- compose_GPGx8cTerm(aa_pool_x8))
parentIonMass(smp.peptide)

pim.GPGx8cTerm <- unlist(lapply(peptides.GPGx8cTerm, function(x){parentIonMass(x)}))
pim.GPx10R <- unlist(lapply(peptides.GPx10R, function(x){parentIonMass(x)}))
pim.iRT <-  unlist(lapply(as.character(iRTpeptides$peptide), function(x){parentIonMass(x)}))

Draw parent ion mass histogram

(pim.min <- min(pim.GPGx8cTerm, pim.GPx10R))
(pim.max <- max(pim.GPGx8cTerm, pim.GPx10R))

(pim.breaks <- seq(round(pim.min - 1) , round(pim.max + 1) , length=75))

hist(pim.GPGx8cTerm, breaks=pim.breaks, probability = TRUE, 
     col='#1111AAAA', xlab='peptide mass [Dalton]', ylim=c(0, 0.006))
hist(pim.GPx10R, breaks=pim.breaks,
     probability = TRUE, add=TRUE, col='#11AA1188')
abline(v=pim.iRT, col='grey')
legend("topleft", c('GPGx8cTerm', 'GPx10R', 'iRT'), 
     fill=c('#1111AAAA', '#11AA1133', 'grey'))

Hydrophobicity

Compute Hydrophobicity value using SSRC

the SSRC model, see @pmid15238601, is implemented as ssrc function in r CRANpkg("protViz").

For a sanity check we apply the ssrc function to a real world LC-MS run peptideStd consits of a digest of the FETUIN_BOVINE protein (400 amol) shipped with r Biocpkg("specL") @pmid25712692.

library(specL)
ssrc <- sapply(peptideStd, function(x){ssrc(x$peptideSequence)})
rt <- unlist(lapply(peptideStd, function(x){x$rt}))
plot(ssrc, rt); abline(ssrc.lm <- lm(rt ~ ssrc), col='red'); 
legend("topleft", paste("spearman", round(cor(ssrc, rt, method='spearman'),2)))

here we apply ssrc to the simulated flycodes and iRT peptides @pmid22577012.

hyd.GPGx8cTerm <- ssrc(peptides.GPGx8cTerm)
hyd.GPx10R <- ssrc(peptides.GPx10R)
hyd.iRT <- ssrc(as.character(iRTpeptides$peptide))

(hyd.min <- min(hyd.GPGx8cTerm, hyd.GPx10R))
(hyd.max <- max(hyd.GPGx8cTerm, hyd.GPx10R))
hyd.breaks <- seq(round(hyd.min - 1) , round(hyd.max + 1) , length=75)

Draw hydrophobicity histogram

hist(hyd.GPGx8cTerm, breaks = hyd.breaks, probability = TRUE, 
     col='#1111AAAA', xlab='hydrophobicity', 
     ylim=c(0, 0.06),
     main='Histogram')
hist(hyd.GPx10R, breaks = hyd.breaks, probability = TRUE, add=TRUE, col='#11AA1188')
  abline(v=hyd.iRT, col='grey')
legend("topleft", c('GPGx8cTerm', 'GPx10R', 'iRT'),  fill=c('#1111AAAA', '#11AA1133', 'grey'))

Quality Control (QC)

QC of composed peptides

Input

round(table(aa_pool_x8)/length(aa_pool_x8), 2)

Output

peptide2aa <- function(seq, from=4, to=4+8){
  unlist(lapply(seq, function(x){strsplit(substr(x, from, to), '')}))
}
peptides.GPGx8cTerm.aa <- peptide2aa(peptides.GPGx8cTerm)
round(table(peptides.GPGx8cTerm.aa)/length(peptides.GPGx8cTerm.aa), 2)
peptides.GPx10R.aa <- peptide2aa(peptides.GPx10R, from=3, to=12)
round(table(peptides.GPx10R.aa)/length(peptides.GPx10R.aa), 2)

Count GP patterns

sample.size 

length(grep('^GP(.*)GP(.*)R$', peptides.GPGx8cTerm))

length(grep('^GP(.*)GP(.*)R$', peptides.GPx10R))

Compute AA frequency table

count the peptides having the same AA composition

sample.size 

table(table(tt<-unlist(lapply(peptides.GPGx8cTerm, 
  function(x){paste(sort(unlist(strsplit(x, ''))), collapse='')}))))
# write.table(tt, file='GPGx8cTerm.txt')
table(table(unlist(lapply(peptides.GPx10R, 
  function(x){paste(sort(unlist(strsplit(x, ''))), collapse='')}))))

the r Biocpkg("NestLink") function plot_in_silico_LCMS_map graphs the LC-MS maps.

par(mfrow=c(2, 2))
h <- NestLink:::.plot_in_silico_LCMS_map(peptides.GPGx8cTerm, main='GPGx8cTerm')
h <- NestLink:::.plot_in_silico_LCMS_map(peptides.GPx10R, main='GPx10R')

Session info

Here is the output of sessionInfo() on the system on which this document was compiled:

sessionInfo()

References



cpanse/NestLink documentation built on May 16, 2022, 2:33 a.m.