Hey there,
hopefully this question is placed correctly on this board. I am programming a targeted analysis for 13C-labeled metabolites using xcms. For the calculation of the accurate masses to look for I wanted to use Rdisop. However, since this library considers only natural isotope distributions, the accurate masses of labeled compounds occurring naturally at very low abundances are not calculated. A short example calculating the accurate mass of a fully 13C-labeled Glucose demonstrates this:
> library( "Rdisop" )
> getIsotope( getMolecule( "C6H12O6" ), 7 )
[1] 186.0774 2.678176e-07
H1 <- 1.007825032
C13 <- 13.00335484
O16 <- 15.99491462
6*C13 + 12 * H1 + 6 * O16
[1] 186.0835
Obviously, the mass should be 186.0835u and not 186.0774u (difference around 30ppm). I assume, that the calculated mass results from a mix of naturally occurring 12C/13C + 16O/18O since the calculated abundance (2.7e-07) is also magnitudes larger than the one of a naturally occurring [13C]6[1H]12[16O]6 which should be 6.4e-10.
My question is if anyone can give me a hint how to calculate the accurate masses for isotopic labeled compounds using Rdisop? The isotope distribution seams to be hard coded in initializePSE(). Is there any hook, to which I could provide an alternative isotope distribution?
Many thanks in advance
Isam
Try:
?initializeElements
You should be able to define your own elements here.
Try initializeCHNOPSMgKCaFe() to see how the list of elements should be build.
Hope it helps.
Hey,
thanks for your reply.
Thats how I started. However, initializeElements just calls Rdisop:::.getElement which again calls initializePSE. From my understanding the initialize methods only allow for the selection of relevant elements but not for the adjustment of isotope distributions.
EDIT:
OK sorry, this was the pointer I needed. I can write my own initialize13C() method, let it call initializePSE, adjust the 13C distribution.
getMolecule("C6H12O6", elements = initialize13C() )
should finally do the trick.
Many thanks!
Isam
initializeCHNOPSMgKCaFe() was just to show you how to build the list of elements you should use with initializeElements.
EDIT: dang... too slow... and I guess initializeElements is not needed if you use my example.
Just for the completeness, this solution works pretty well:
initialize13C <- function( a12C13C14C = c( 0.1, 0.9, 0 ) ) {
require("Rdisop")
elements <- initializePSE()
elements <- lapply( elements, function(x){
if(x$name == "C" ){
x$isotope$abundance <- a12C13C14C
}
return( x )
})
return(elements)
}
H1 <- 1.007825032
C12 <– 12
C13 <- 13.00335484
O16 <- 15.99491462
all.equal(
sapply( 0:6, FUN=function(x) { x * C13 + (6-x) * C12 + 12 * H1 + 6 * O16 } ),
getIsotope( getMolecule( "C6H12O6", elements=initialize13C() ) )[ 1,1:7 ]
)
[1] TRUE
Thanks again for pointing me to the right direction. I simply missed the "elements" parameter in getMolecule().
Isam