4 Workflow
(Note: The following examples are for illustrative purpose only and need not be biologically meaningful.)
Once both the tumor genome(s) and the given mutational signatures have been read (see above), the contribution of the given signatures to the somatic mutations in individual tumors can be determined using the following workflow.
In the following, we assume to have the signatures in a list object named signatures and the mutation frequencies of the tumor genome(s) in another list object named genomes.
Important note: it is imperative that the mutation frequencies represented by both signatures and genomes have been computed with the same characteristics. That is, if the signatures refer to sequences of 5 nucleotides/bases (with the mutated base at the center), then also the genomes must have been read for 5 bases. If the signatures have been produced taking the transcription direction into account, then also the genomes must have been read taking the transcription direction into account, and so on.
4.1 Visualizing genome characteristics and mutational signatures
Given that the signatures and the genomes (if read appropriately) have the same format and contain the same type of information (fractions/percentages of somatic mutations that have specific features, e.g., specific neighboring bases), both can essentially be visualized in the same way.
The function plotMutationDistribution() takes as an input either a single signature or the mutation frequencies of an individual tumor genome (either of Alexandrov- or of Shiraishi-type) and plots the mutation frequency data according to the signature model.
To plot, for example, Alexandrov/COSMIC signature 1 (obtained as described in Section 3.1.1):
signatures <- readAlexandrovSignatures()## Warning in file(file, "rt"): "internal" method cannot handle https
## redirection to: 'https://cancer.sanger.ac.uk/cancergenome/assets/
## signatures_probabilities.txt'## Warning in file(file, "rt"): "internal" method failed, so trying "libcurl"plotMutationDistribution(signatures[[1]])
To plot one of the Shiraishi signatures provided with this package (see Section 3.1.2):
sigfiles <- system.file("extdata",
paste0("Nik-Zainal_PMID_22608084-pmsignature-sig",1:4,".tsv"),
package="decompTumor2Sig")
signatures <- readShiraishiSignatures(files=sigfiles)
plotMutationDistribution(signatures[[1]])
decompTumor2Sig’s representation of the mutation frequency data of Shiraishi-type signatures uses sequence logos for the flanking bases and the variant bases (with the heights of the bases being proportional to their probability/frequency). The two possibilities for the mutated central base (C or T) are represented next to each other and their respective frequency is indicated below. This side-by-side representation allows to distinguish the probabilities of variant bases (on top) according to the mutated base. Transcription strand bias (if information on transcription direction is used) is shown in the upper right corner (frequency of mutations on the transcription strand, “+”, and the opposite strand, “-”).
In the plot above, for example, C and T are nearly equally likely to be mutated by the represented mutational process, but a mutated C most frequently becomes a T, while a mutated T becomes one of the other bases with roughly equal probability. Also, the mutational signature has next to no strand bias.
This representation is similar to the way the pmsignature package represents such signatures, as shown by the following example:
(This plot above was generated with pmsignature and serves only for comparison, showing the same signature as above.)
To show that genome mutation frequencies can be represented in the same manner, the following example reads the tumor genomes provided with this package using the Alexandrov model (see Section 3.2.1) and plots the mutation frequencies of the first genome:
refGenome <- BSgenome.Hsapiens.UCSC.hg19::BSgenome.Hsapiens.UCSC.hg19
gfile <- system.file("extdata",
"Nik-Zainal_PMID_22608084-VCF-convertedfromMPF.vcf.gz",
package="decompTumor2Sig")
genomes <- readGenomesFromVCF(gfile, numBases=3, type="Alexandrov",
trDir=FALSE, refGenome=refGenome, verbose=FALSE)
plotMutationDistribution(genomes[[1]])
4.2 Explained variance as a function of the number of signatures
In many cases already a small subset of the given signatures is sufficient to explain the major part of the variance of the mutation frequencies observed in a single tumor genome.
The explained variance can be determined by comparing the true mutation frequencies \(g_i\) of the tumor genome—where \(i\) is one of the mutation features of the signature model (e.g., triplet mutations for Alexandrov signatures, or base changes or flanking bases for Shiraishi signatures)—to the mutation frequencies \(\hat{g}_i\) obtained when re-composing/reconstructing the mutation frequencies of the tumor genome from the mutational signatures and their computed exposures/contributions. (See Section 4.3 for computing the exposures/contributions, and Section 4.4 for reconstructing the mutation frequencies of a tumor.)
For the Alexandrov model, the explained variance can be estimated as:
\[ \mathrm{evar}=1-\frac{\sum_i{(g_i-\hat{g}_i)^2}}{\sum_i({g_i}-\bar{g})^2} \]
where the numerator is the residual sum of squares (RSS) between the predicted and true mutation frequencies of the tumor genome (i.e., the squared error), and the denominator can be interpreted as the deviation from a tumor genome with a uniform mutation frequency of 1/96 for each feature (which is precisely the average mutation frequency \(\bar{g}\)).
For the Shiraishi model, \(\bar{g}\) does not describe a tumor genome with uniform mutation frequencies, and hence the explained variance is estimated as:
\[\mathrm{evar}=1-\frac{\sum_i{(g_i-\hat{g}_i)^2}}{\sum_i({g_i}-{g}^{*}_i)^2}\]
where \(g^{*}\) is a uniform tumor model that uses mutation frequencies of 1/6 for the six possible base changes, 1/4 for each of the possible flanking bases, and 1/2 for each of the two possible transcription-strand directions. For more details, please see Krüger and Piro (2018).
The function plotExplainedVariance() allows to visually analyze how many signatures are necessary to explain certain fractions of the variance of a tumor genome’s mutational patterns.
For an increasing number K of signatures, the highest variance explained by subsets of K signatures will be plotted. This can help to evaluate what minimum threshold for the explained variance could be used to decompose tumor genomes with the function decomposeTumorGenomes() (see below).
4.2.1 Example: input data
As a simple example, load a set of 15 Shiraishi-type signatures (object signatures) provided with this package. These signatures were obtained with the package pmsignature from a set of 435 tumor genomes with at least 100 somatic mutations from ten different tumor entities (data from Alexandrov et al. Nature 500:415-421, 2013; for the analysis, see our papers mentioned in Section 1.1):
# load the 15 Shiraishi signatures obtained from
# 435 tumor genomes from Alexandrov et al.
sfile <- system.file("extdata",
"Alexandrov_PMID_23945592_435_tumors-pmsignature-15sig.Rdata",
package="decompTumor2Sig")
load(sfile)
length(signatures)## [1] 15signatures[1]## [[1]]
## [,1] [,2] [,3] [,4] [,5] [,6]
## [1,] 6.431522e-13 2.320704e-02 2.293643e-79 0.64050260 0.1876257 0.1486647
## [2,] 2.916982e-01 2.059223e-01 2.043698e-01 0.29800971 0.0000000 0.0000000
## [3,] 2.239856e-01 4.126288e-01 1.982799e-01 0.16510571 0.0000000 0.0000000
## [4,] 9.822384e-02 7.036368e-14 8.147461e-01 0.08703009 0.0000000 0.0000000
## [5,] 2.524624e-01 1.965627e-01 2.681366e-01 0.28283839 0.0000000 0.0000000
## [6,] 4.422814e-01 5.577186e-01 0.000000e+00 0.00000000 0.0000000 0.0000000This loads an object signatures with 15 Shiraishi signatures for mutated subsequences of length 5 (five nucleotides with the mutated base at the center) and taking transcription direction into account.
Now read the tumor genomes (object genomes) provided with this package, as described in Section 3.2.1:
# load the reference genome and the transcript annotation database
refGenome <- BSgenome.Hsapiens.UCSC.hg19::BSgenome.Hsapiens.UCSC.hg19
transcriptAnno <-
TxDb.Hsapiens.UCSC.hg19.knownGene::TxDb.Hsapiens.UCSC.hg19.knownGene
# take the breast cancer genomes from Nik-Zainal et al (PMID: 22608084)
gfile <- system.file("extdata",
"Nik-Zainal_PMID_22608084-VCF-convertedfromMPF.vcf.gz",
package="decompTumor2Sig")
# read the cancer genomes in VCF format
genomes <- readGenomesFromVCF(gfile, numBases=5, type="Shiraishi",
trDir=TRUE, refGenome=refGenome,
transcriptAnno=transcriptAnno, verbose=FALSE)4.2.2 Example: plot the explained variance
The explained variance can be plotted only for a single tumor genome. Plotting the explained variance of the first tumor genome when using increasing subsets of the 15 mutational signatures (see above), for example, can be done with the following command:
plotExplainedVariance(genomes[[1]], signatures, minExplainedVariance=0.9,
minNumSignatures=2, maxNumSignatures=NULL,
greedySearch=TRUE)
The function plotExplainedVariance() takes the following arguments:
- genome: The mutation frequencies of a single tumor genome (according to the Alexandrov or Shiraishi model).
- signatures: The set of given signatures (must be of the same model as the genome).
- minExplainedVariance (default is NULL): If specified, the smallest subset of signatures necessary to explain at least this fraction of the variance is highlighted in red (including the list of the signatures in the subset).
- minNumSignatures (default is 2): Take at least this minimum number of signatures.
- maxNumSignatures (default is NULL, i.e., all): If specified, take at most this maximum number of signatures.
- greedySearch (default is FALSE): If FALSE, the function will evaluate for each number of signatures K all possible subsets of K signatures to compute the highest explained variance for K. This can take very long if the total number of given signatures is too high. If TRUE, first, all possible subsets with K=minNumSignatures are evaluated, taking the subset with the highest explained variance, then stepwise one additional signature (highest increase in explained variance) is added to the already identified set. This approximation is much faster but not guaranteed to plot the maximum explained variance for all K (a greedy search can get stuck in a local optimum).
4.3 Decomposing tumor genomes by signature refitting (contribution prediction)
This is the heart of the functionality of decompTumor2Sig, its main purpose.
Given a tumor genome and a set of mutational signatures (that represent mutational processes like UV light, smoking, etc; see Alexandrov and Stratton, Curr Opin Genet Dev 24:52-60, 2014), we would like to estimate how strongly the different signatures (processes) contributed to the overall mutational load observed in the tumor. To this end decompTumor2Sig relies on quadratic programming. For details, see our paper (Krüger and Piro, 2018).
The result will be a vector of weights/contributions (or “exposures”) which indicate the fractions or percentages of somatic mutations which can likely be attributed to the given signatures.
Lets take, for instance, the signature and tumor data used for the example in Section 4.2 (see there).
To compute the contributions of the 15 given signatures to the first tumor genome, we can run the following command …
exposure <- decomposeTumorGenomes(genomes[1], signatures)… and get the following exposures (contributions) for the 15 signatures:
exposure## $PD3851a
## sign_1 sign_2 sign_3 sign_4 sign_5
## 4.874215e-02 1.237708e-01 5.984978e-02 5.764690e-02 1.426703e-20
## sign_6 sign_7 sign_8 sign_9 sign_10
## 1.109195e-01 8.766150e-02 5.115825e-02 2.758243e-02 9.180615e-02
## sign_11 sign_12 sign_13 sign_14 sign_15
## 3.070582e-02 3.256016e-02 3.608152e-02 1.953982e-01 4.611685e-02The exposures/contributions for a single tumor genome can also be plotted:
plotDecomposedContribution(exposure)
In some cases, multiple tumor genomes need to be decomposed (each one, however, individually). In this case, a set of tumor genomes can be passed to decomposeTumorGenomes():
exposures <- decomposeTumorGenomes(genomes, signatures)
length(exposures)## [1] 21exposures[1:2]## $PD3851a
## sign_1 sign_2 sign_3 sign_4 sign_5
## 4.874215e-02 1.237708e-01 5.984978e-02 5.764690e-02 1.426703e-20
## sign_6 sign_7 sign_8 sign_9 sign_10
## 1.109195e-01 8.766150e-02 5.115825e-02 2.758243e-02 9.180615e-02
## sign_11 sign_12 sign_13 sign_14 sign_15
## 3.070582e-02 3.256016e-02 3.608152e-02 1.953982e-01 4.611685e-02
##
## $PD3890a
## sign_1 sign_2 sign_3 sign_4 sign_5 sign_6
## 0.01525397 0.02791685 0.07609293 0.07239998 0.00815879 0.17830407
## sign_7 sign_8 sign_9 sign_10 sign_11 sign_12
## 0.04543191 0.09344739 0.14291157 0.13844086 0.05390118 0.00000000
## sign_13 sign_14 sign_15
## 0.03827808 0.08551768 0.02394473decompTumor2Sig provides an additional function that can be used to verify whether a list object is of the same format as it is returned by decomposeTumorGenomes():
isExposureSet(exposures)## [1] TRUE4.3.1 Finding a subset of signatures with a minimum explained variance
Instead of decomposing a tumor genome precisely into the given set of signatures (15 in the example above), the function decomposeTumorGenomes() can alternatively be used to search for subsets of signatures for which the decomposition satisfies a minimum threshold of explained variance (here, we show this only for the first tumor):
exposures <- decomposeTumorGenomes(genomes[1], signatures,
minExplainedVariance=0.9,
minNumSignatures=2, maxNumSignatures=NULL,
greedySearch=FALSE, verbose=TRUE)## Decomposing genome 1 (PD3851a) with 2 signatures ...
## with 3 signatures ...
## with 4 signatures ...
## with 5 signatures ...
## with 6 signatures ...
## with 7 signatures ...
## explained variance: 0.903349612580512For each tumor genome, the minimum subset of signatures that explain at least minExplainedVariance percent of the variance of the mutation frequencies will be identified (the exposures of all other signatures will be set to NA):
exposures## $PD3851a
## sign_1 sign_2 sign_3 sign_4 sign_5 sign_6
## NA 0.14986576 0.08839117 0.11203026 NA 0.18900885
## sign_7 sign_8 sign_9 sign_10 sign_11 sign_12
## 0.16543502 NA NA 0.09461580 NA NA
## sign_13 sign_14 sign_15
## NA 0.20065314 NAplotDecomposedContribution(exposures[[1]])
[If plotDecomposedContribution() is run with removeNA=FALSE, also signatures with an NA value as exposure will be included in the x-axis of the plot. Additionally, the signatures can be passed to the function using the parameter signatures; if so, signature names for the plot will be taken from this object, otherwise they are inferred from the exposure object or set sign_1 to sign_N.]
Important note: although a subset of signatures may explain the somatic mutations observed in the tumor genomes reasonably well, they need not necessarily be the signatures with the highest contribution when taking the entire set of signatures (e.g., due to a greedy search which can get stuck in a local optimum).
Important note: if for a tumor genome no (sub)set of signatures is sufficient to explain at least minExplainedVariance percent of the variance, no result (NULL) is returned for that tumor.
The search behavior of decomposeTumorGenomes() when finding a subset of signatures to explain the somatic mutations of a tumor can be influenced by the following arguments:
- minExplainedVariance (default is NULL): If not specified, exactly maxNumSignatures (see below; default: all) will be taken for decomposing each genome; if specified (between 0 and 1), the smallest subset of signatures which explains at least minExplainedVariance of the variance is taken for the decomposition.
- minNumSignatures (default is 2): If a search for a subset of signatures is performed, at least minNumSignatures will be taken.
- maxNumSignatures (default is NULL, i.e., all): If a search for a subset of signatures is performed, at most maxNumSignatures will be taken; if NULL, all given signatures will be taken as the maximum; if maxNumSignatures is specified but minExplainedVariance=NULL (no search), then exactly the best maxNumSignatures will be taken.
- greedySearch (default is FALSE): If FALSE, a full search will be performed: for increasing numbers K of signatures, all possible subsets of K signatures will be tested and the subset with the highest explained variance is chosen, increasing K until the threshold of minExplainedVariance is satisfied; if TRUE, a much faster, greedy search is performed: first, all possible subsets with K=minNumSignatures are evaluated, taking the subset with the highest explained variance. Then, stepwise one additional signature at a time (highest increase in explained variance) is added to the already identified set until the threshold of minExplainedVariance is satisfied.
Performing, for example, a greedy search for the example above is much faster (shown only for the first tumor):
exposures <- decomposeTumorGenomes(genomes[1], signatures,
minExplainedVariance=0.95,
minNumSignatures=2, maxNumSignatures=NULL,
greedySearch=TRUE, verbose=TRUE)## Decomposing genome 1 (PD3851a) with 2 signatures ...
## adding signature 3 ...
## adding signature 4 ...
## adding signature 5 ...
## adding signature 6 ...
## adding signature 7 ...
## adding signature 8 ...
## adding signature 9 ...
## adding signature 10 ...
## explained variance: 0.959544941532486exposures## $PD3851a
## sign_1 sign_2 sign_3 sign_4 sign_5 sign_6
## 0.04423271 0.17855158 0.05436733 0.07616399 NA 0.12658971
## sign_7 sign_8 sign_9 sign_10 sign_11 sign_12
## 0.17677283 0.03854308 0.03481150 0.07507863 NA NA
## sign_13 sign_14 sign_15
## NA 0.19488865 NAplotDecomposedContribution(exposures[[1]])
Important note: of course, a greedy search which starts from the best combination of minNumSignatures need not yield the same result as a full search because the latter finds the overall best subset while the greedy search can get stuck in a local optimum depending on the starting point of the search. The precise behavior depends on different factors, including for example the similarity between mutational signatures and the minimum required explained variance (lower thresholds can easily be satisfied by very different combinations of signatures). We recommend to test different settings (minimum numbers of signatures, thresholds for explained variance, etc) and learn about the behavior to ensure a meaningful biological interpretation of the results.
4.3.2 Computing the explained variance
Given the mutation frequencies of one or more tumor genomes (genomes), a set of mutational signatures (signatures) and their computed exposures/contributions to the given tumor (exposures), the following command allows to compute—for each individual tumor—the variance of the mutation frequency data that the exposures explain.
The following is an example taking the full set of tumors and signatures from Section 4.2.1:
exposures <- decomposeTumorGenomes(genomes, signatures)
computeExplainedVariance(exposures, signatures, genomes)## PD3851a PD3890a PD3904a PD3905a PD3945a PD4005a PD4006a
## 0.9796806 0.9693211 0.9913913 0.9695197 0.9861931 0.9956940 0.8988394
## PD4085a PD4086a PD4088a PD4103a PD4107a PD4109a PD4115a
## 0.9430757 0.9912875 0.9693047 0.9843990 0.9910398 0.9931142 0.9607950
## PD4116a PD4120a PD4192a PD4194a PD4198a PD4199a PD4248a
## 0.9937769 0.9995375 0.9764347 0.9678732 0.9866895 0.9987643 0.9841221Note: for computing the variance explained for a single tumor by a set of signatures, the corresponding number of exposure values must be the same as the number of signatures. Undefined exposure values (NA), which can be present if a search for a subset of signatures has been performed as described above, will be set to zero such that the corresponding signature does not contribute. For genomes for which the minExplainedVariance could not be reached, and whose exposure vectors are NULL, the explained variance will be set to NA.
4.4 Re-composing/reconstructing tumor genomes from exposures and signatures
Estimating the explained variance of the decomposition of a tumor genome (see Sections 4.2 and 4.3.2) and assessing its quality requires the mutation frequencies \(\hat{g}_i\) of the tumor genome to be reconstructed, or predicted, from the mutational signatures \(S_j\) and their exposures/contributions (or weights) \(w_j\):
\[ \hat{g}_i = \sum_j{w_j * (S_j)_i} \]
This can be easily achieved using the function composeGenomesFromExposures(). The following is an example taking the tumor genomes from Section 3.2, the signatures from Section 4.2.1, and the exposures as computed in Section 4.3:
genomes_predicted <- composeGenomesFromExposures(exposures, signatures)
genomes_predicted[1:2]## $PD3851a
## [C>A] [C>G] [C>T] [T>A] [T>C] [T>G]
## mut 0.1871191 0.1003964 0.3534400 0.1183363 0.1637962 0.07691213
## -2 0.2789142 0.2232669 0.1935118 0.3043071 0.0000000 0.00000000
## -1 0.2759081 0.2481861 0.1757908 0.3001150 0.0000000 0.00000000
## +1 0.2757096 0.1900909 0.2655066 0.2686929 0.0000000 0.00000000
## +2 0.2825315 0.2121257 0.2028239 0.3025190 0.0000000 0.00000000
## tr 0.4704991 0.5295009 0.0000000 0.0000000 0.0000000 0.00000000
##
## $PD3890a
## [C>A] [C>G] [C>T] [T>A] [T>C] [T>G]
## mut 0.1647946 0.2312420 0.2308237 0.1193547 0.136055 0.1177298
## -2 0.2896184 0.2248804 0.1791938 0.3063074 0.000000 0.0000000
## -1 0.2353792 0.2096112 0.1716452 0.3833645 0.000000 0.0000000
## +1 0.2882678 0.2100574 0.1572802 0.3443946 0.000000 0.0000000
## +2 0.3051997 0.1888451 0.1982830 0.3076722 0.000000 0.0000000
## tr 0.4883765 0.5116235 0.0000000 0.0000000 0.000000 0.0000000Once genomes have been reconstructed, the function evaluateDecompositionQuality() allows to compare the reconstructed genome mutation features to the originally observed features in order to assess the quality of the decomposition. The function is applied to an individual tumor genome, and can either return numerical quality measurements, or provide a quality plot which includes said measurements.
Numerical measurements to compare the reconstructed/predicted and the original tumor mutation patterns are:
- explainedVariance: The explained variance (see Section 4.2).
- pearsonCorr: The Pearson correlation coefficient (PCC) between the predicted/reconstructed and the original mutation frequencies of the tumor genome. Although the PCC does usually not consider the amplitudes of two input data vectors—only their linear relationship—here, a high PCC also entails small absolute differences in the mutation frequencies because the frequencies in both the predicted and the original data sum up to 1 (i.e., they are normalized at the same level), such that a high correlation automatically means very similar values.
Example for obtaining numerical measurements:
evaluateDecompositionQuality(exposures[[1]], signatures,
genomes[[1]], plot=FALSE)## $explainedVariance
## [1] 0.9796806
##
## $pearsonCorr
## [1] 0.9989003Example for obtaining a quality plot which compares the reconstructed and the original data:
evaluateDecompositionQuality(exposures[[1]], signatures,
genomes[[1]], plot=TRUE)
4.5 Mapping and comparing sets of signatures
In some cases it may be of interest to compare or find a mapping between two sets of signatures, e.g., if they have been inferred from different datasets. For this purpose, decompTumor2Sig provides a set of additional functions.
4.5.1 Comparison of signatures of the same format
Let’s first read two distinct sets of Shiraishi signatures of the same format (5 bases, with transcript-strand direction):
# get 4 Shiraishi signatures from 21 breast cancers from
# Nik-Zainal et al (PMID: 22608084)
sigfiles <- system.file("extdata",
paste0("Nik-Zainal_PMID_22608084-pmsignature-sig",1:4,".tsv"),
package="decompTumor2Sig")
sign_s4 <- readShiraishiSignatures(files=sigfiles)
# get 15 Shiraishi signatures obtained from
# 435 tumor genomes from Alexandrov et al.
sfile <- system.file("extdata",
"Alexandrov_PMID_23945592_435_tumors-pmsignature-15sig.Rdata",
package="decompTumor2Sig")
load(sfile)
sign_s15 <- signaturesSince these signatures have the same format, we can directly compare them. Using the function determineSignatureDistances(), we can compute the distances of one given signature to all target signatures from a set:
determineSignatureDistances(fromSignature=sign_s4[[1]], toSignatures=sign_s15,
method="frobenius")## sign_1 sign_2 sign_3 sign_4 sign_5 sign_6 sign_7
## 0.7585028 1.1070702 1.1324959 0.7741531 1.2625974 0.9338330 0.6743049
## sign_8 sign_9 sign_10 sign_11 sign_12 sign_13 sign_14
## 1.0077730 1.5575345 1.2758753 1.1879876 1.2602450 1.0191394 0.6792248
## sign_15
## 1.0914233Apart from the Frobenius distance (method=“frobenius”), which is suitable to compare matrices and hence Shiraishi signatures, other distance metrics can be used: the “rss” (residual sum of squares = squared error) or any distance measure available for the function dist of the package stats.
If not the distances of one signature to an entire signature set is needed, but instead a mapping from one signature set to another, mapSignatureSets() can be used.
mapSignatureSets(fromSignatures=sign_s4, toSignatures=sign_s15,
method="frobenius", unique=FALSE)## Nik-Zainal_PMID_22608084-pmsignature-sig1.tsv
## "sign_7"
## Nik-Zainal_PMID_22608084-pmsignature-sig2.tsv
## "sign_9"
## Nik-Zainal_PMID_22608084-pmsignature-sig3.tsv
## "sign_6"
## Nik-Zainal_PMID_22608084-pmsignature-sig4.tsv
## "sign_12"Like for determineSignatureDistances(), with the function mapSignatureSets() a mapping can be built based on different distance metrics. Additionally, the user can specify whether the mapping should be unique (one-to-one mapping), or not.
If unique=FALSE then for each signature of fromSignatures the best match (minimum distance) of toSignatures is selected. The selected signatures need not be unique, i.e., one signature of toSignatures may be the best match for multiple signatures of fromSignatures.
If unique=TRUE, i.e., if a unique (one-to-one) mapping is required, an iterative procedure is performed: in each step, the best matching pair from fromSignatures and toSignatures is mapped and then removed from the list of signatures that remain to be mapped, such that they cannot be selected again. In this case, of course, fromSignatures must not contain more signatures than toSignatures.
4.5.2 Comparison of signatures of different types or formats
Sometimes it may also be useful to compare different types or formats of signatures. For example, since Alexandrov signatures are comparably well studied, it might be interesting to determine which Alexandrov signature is closest to a Shiraishi signature of interest.
Since only signatures of the same format can be directly compared or mapped (see above), decompTumor2Sig provides two functions that transform signatures, such that two signatures, or two sets of signatures, can be converted to the same format.
One of these functions, convertAlexandrov2Shiraishi(), has already been presented in Section 3.1.4 (see there for details). We can convert Alexandrov signatures to Shiraishi-type signatures with 3 bases (without transcription-strand direction):
sign_a <- readAlexandrovSignatures()
sign_a2s <- convertAlexandrov2Shiraishi(sign_a)Note: Of course, there is some information loss here (better: a loss of specificity), as we discuss in our paper (Krüger and Piro, 2018).
Additionally, the function downgradeShiraishiSignatures() can be used to reduce the number of flanking bases and/or discard the information on transcription-strand direction from one or more Shiraishi signatures:
sign_sdown <- downgradeShiraishiSignatures(sign_s15, numBases=3,
removeTrDir=TRUE)
sign_s15[1]## [[1]]
## [,1] [,2] [,3] [,4] [,5] [,6]
## [1,] 6.431522e-13 2.320704e-02 2.293643e-79 0.64050260 0.1876257 0.1486647
## [2,] 2.916982e-01 2.059223e-01 2.043698e-01 0.29800971 0.0000000 0.0000000
## [3,] 2.239856e-01 4.126288e-01 1.982799e-01 0.16510571 0.0000000 0.0000000
## [4,] 9.822384e-02 7.036368e-14 8.147461e-01 0.08703009 0.0000000 0.0000000
## [5,] 2.524624e-01 1.965627e-01 2.681366e-01 0.28283839 0.0000000 0.0000000
## [6,] 4.422814e-01 5.577186e-01 0.000000e+00 0.00000000 0.0000000 0.0000000sign_sdown[1]## [[1]]
## [C>A] [C>G] [C>T] [T>A] [T>C] [T>G]
## mut 6.431522e-13 2.320704e-02 2.293643e-79 0.64050260 0.1876257 0.1486647
## -1 2.239856e-01 4.126288e-01 1.982799e-01 0.16510571 0.0000000 0.0000000
## +1 9.822384e-02 7.036368e-14 8.147461e-01 0.08703009 0.0000000 0.0000000Having obtained two signature sets of the same format (sequence triplets, but no transcription-strand direction), we can now map one set to the other:
mapSignatureSets(fromSignatures=sign_sdown, toSignatures=sign_a2s,
method="frobenius", unique=TRUE)## sign_1 sign_2 sign_3 sign_4 sign_5
## "Signature.22" "Signature.4" "Signature.10" "Signature.5" "Signature.15"
## sign_6 sign_7 sign_8 sign_9 sign_10
## "Signature.8" "Signature.26" "Signature.28" "Signature.13" "Signature.30"
## sign_11 sign_12 sign_13 sign_14 sign_15
## "Signature.24" "Signature.2" "Signature.18" "Signature.1" "Signature.12"