peptideLevel_PresAbsDE: Presence/Absence peptide-level analysis

In YuliyaLab/ProteoMM: Multi-Dataset Model-based Differential Expression Proteomics Analysis Platform

Presence/Absence peptide-level analysis

Description

Presence/Absence peptide-level analysis uses all peptides for a protein as IID to produce 1 p-value across multiple (2+) datasets. Significance is estimated using a g-test which is suitable for two treatment groups only.

Usage

peptideLevel_PresAbsDE(mm, treatment, prot.info, pr_ppos = 2)

Arguments

mm	m x n matrix of intensities, num peptides x num samples
treatment	vector indicating the treatment group of each sample ie [1 1 1 1 2 2 2 2...]
prot.info	2+ colum data frame of peptide ID, protein ID, etc columns
pr_ppos	- column index for protein ID in prot.info. Can restrict to be #2...

Value

A list of length two items:

ProtIDused

protein identification information taken from prot.info, a column used to identify proteins

FC

Approximation of the fold change computed as percent missing observations group 1 munis in percent missing observations group 2

P_val

p-value for the comparison between 2 groups (2 groups only here)

BH_P_val

Benjamini-Hochberg adjusted p-values

statistic

statistic returned by the g-test, not very useful as depends on the direction of the test and can produce all 0's

num_peptides

number of peptides within a protein

metadata

all columns of metadata from the passed in matrix

Examples

# Load mouse dataset
data(mm_peptides)
head(mm_peptides)
intsCols = 8:13 # different from parameter names as R uses
                # outer name spaces if variable is undefined
metaCols = 1:7 # reusing this variable
m_logInts = make_intencities(mm_peptides, intsCols)  # will reuse the name
m_prot.info = make_meta(mm_peptides, metaCols)
m_logInts = convert_log2(m_logInts)
grps = as.factor(c('CG','CG','CG', 'mCG','mCG','mCG'))

set.seed(135)
mm_m_ints_eig1 = eig_norm1(m=m_logInts,treatment=grps,prot.info=m_prot.info)
mm_m_ints_eig1$h.c # check the number of bias trends detected
mm_m_ints_norm = eig_norm2(rv=mm_m_ints_eig1)

# Load human dataset
data(hs_peptides)
head(hs_peptides)
intsCols = 8:13 # different from parameter names as R
                # uses outer name spaces if variable is undefined
metaCols = 1:7 # reusing this variable
m_logInts = make_intencities(hs_peptides, intsCols)  # will reuse the name
m_prot.info = make_meta(hs_peptides, metaCols)
m_logInts = convert_log2(m_logInts)
grps = as.factor(c('CG','CG','CG', 'mCG','mCG','mCG'))

set.seed(137) # different seed for different organism
hs_m_ints_eig1 = eig_norm1(m=m_logInts,treatment=grps,prot.info=m_prot.info)
hs_m_ints_eig1$h.c # check the number of bias trends detected
hs_m_ints_norm = eig_norm2(rv=hs_m_ints_eig1)

# Set up for presence/absence analysis
raw_list = list()
norm_imp_prot.info_list = list()
raw_list[[1]] = mm_m_ints_eig1$m
raw_list[[2]] = hs_m_ints_eig1$m
norm_imp_prot.info_list[[1]] = mm_m_ints_eig1$prot.info
norm_imp_prot.info_list[[2]] = hs_m_ints_eig1$prot.info

protnames_norm_list = list()
protnames_norm_list[[1]] = unique(mm_m_ints_norm$normalized$MatchedID)
protnames_norm_list[[2]] = unique(hs_m_ints_norm$normalized$MatchedID)

presAbs_dd = get_presAbs_prots(mm_list=raw_list,
                              prot.info=norm_imp_prot.info_list,
                              protnames_norm=protnames_norm_list,
                              prot_col_name=2)

presAbs_de = peptideLevel_PresAbsDE(presAbs_dd[[1]][[1]],
                                    grps, presAbs_dd[[2]][[1]],
                                    pr_ppos=2)

YuliyaLab/ProteoMM documentation built on April 19, 2022, 8:12 a.m.