DCCM.XYZ using LMI method

[Former user]

I am trying to perform dynamic cross correlation on 167 residues over 10000 frames.
Script:

trj <- read.dcd("temp.dcd")

pdb <- read.pdb("temp.pdb")

ca.inds <- atom.select(pdb, elety="CA")

xyz <- fit.xyz(fixed=pdb$xyz, mobile=trj, fixed.inds=ca.inds$xyz, mobile.inds=ca.inds$xyz,ncore = 4, nseg.scale = 4)

cij<-dccm(xyz[,ca.inds$xyz], ncore=4, nseg.scale=4) cov2dccm(vcov,method="lmi", ncore=4)

Error

NATOM = 2648 NFRAME= 10001 ISTART= 0 last = 10001 nstep = 10001 nfile = 10001 NSAVE = 1 NDEGF = 0 version 24 |======================================================================| 100% Error in matrix(nrow = l, ncol = l) : invalid 'nrow' value (too large or NA) Calls: cov2dccm -> matrix Execution halted

[Xinqiu Yao]

Hi,

It is possibly just the data is too "big". Could you reduce the trajectory file, e.g. by skipping every 10 frames (so finally you will get 1,000 frames to analyse), and try it again? Let me know if it still gives errors.

[Rupesh Agarwal]

NATOM = 2648 NFRAME= 1001 ISTART= 0 last = 1001 nstep = 1001 nfile = 1001 NSAVE = 1 NDEGF = 0 version 24 |======================================================================| 100% Error in matrix(nrow = l, ncol = l) : invalid 'nrow' value (too large or NA) Calls: cov2dccm -> matrix Execution halted

same for 1000 frames

NATOM = 2648 NFRAME= 101 ISTART= 0 last = 101 nstep = 101 nfile = 101 NSAVE = 1 NDEGF = 0 version 24 |======================================================================| 100% Error in matrix(nrow = l, ncol = l) : invalid 'nrow' value (too large or NA) Calls: cov2dccm -> matrix Execution halted

and 100 frames

[Rupesh Agarwal]

What is the difference between: dccm(xyz) and dccm(xyz[,ca.inds$xyz], ncore=4, nseg.scale=4) cov2dccm(vcov,method="lmi", ncore=4)? Does this dccm(xyz) uses pearson? And if I want to get dynamic cross correlation using LMI, can I use this without dccm(xyz)?

covmat <- cov(xyz[,ca.inds$xyz])

ccmat <- cov2dccm(covmat,method="lmi", ncore = 4)

[Xinqiu Yao]

dccm(xyz) calculate DCCM (i.e. Dynamical Cross-Correlation Matrix) for all the atoms while dccm(xyz[, ca.inds$xyz]) just for C-alpha atoms.

For LMI, you should use the lmi() function (See help(lmi)), although the codes you provide should give the same results...

To developers: We probably should combine lmi() and dccm() and provide a "method" option in dccm() for users to choose either 'Pearson' or 'LMI'. It will be tidier. Thought?

[Barry Grant]

Sounds reasonable Xinqiu. Not sure many use LMI at the moment...

[Rupesh Agarwal]

dccm(xyz[,ca.inds$xyz], ncore=4, nseg.scale=4) cov2dccm(vcov,method="lmi", ncore=4)

covmat <- cov(xyz[,ca.inds$xyz]) ccmat <- cov2dccm(covmat,method="lmi", ncore = 4)

These are same?

[Xinqiu Yao]

No, they are not the same.

dccm(xyz[,ca.inds$xyz], ncore=4, nseg.scale=4) only calculate Pearson correlation.

covmat <- cov(xyz[,ca.inds$xyz])
ccmat <- cov2dccm(covmat,method="lmi", ncore = 4)

Above two commands collectively calculate lmi correlaiton.

As I said, you can simply call lmi(xyz[, ca.inds$xyz], ncore=4) for lmi correlaiton.

[Rupesh Agarwal]

Thanks XIn.

So default is Pearson (if I do not specify anything in dccm). Is that correct?

The problem that you might want to look into for other users is, when you run this it gives the error I mentioned in the initial post cij<-dccm(xyz[,ca.inds$xyz], ncore=4, nseg.scale=4) cov2dccm(vcov,method="lmi", ncore=4).

cov2dccm works fine for same number of frames.

[Xinqiu Yao]

Why do you call cov2dccm(vcov,method="lmi", ncore=4) right after cij<-dccm(xyz[,ca.inds$xyz], ncore=4, nseg.scale=4)? These are two commands!

[Rupesh Agarwal]

http://thegrantlab.org/bio3d/html/dccm.xyz.html

"dccm"(x, reference = NULL, grpby=NULL, ncore=1, nseg.scale=1, ...) cov2dccm(vcov, method = c("pearson", "lmi"), ncore = NULL)

That's how it shows in the tutorial.

[Rupesh Agarwal]

The default for dccm() is Pearson (if I do not specify anything in dccm). Is that correct?

[Xinqiu Yao]

dccm() only calculates Pearson and currently it does not support other methods (for LMI use lmi()).

[Rupesh Agarwal]

Thanks for the clarification.

[Rupesh Agarwal]

What is the difference between LMI (Linear mutual information) and MI (Mutual information)?

Please explain as it was not clear from the code. thanks

[Xinqiu Yao]

Well, it's hard to explain with just a few sentences. Please read the paper: Lange, O.F. and Grubmuller, H. (2006) Proteins. 62:1053-1061.

Simply, LMI calculates only linear correlation (like Pearson) whereas MI captures also non-linear correlations. Note that bio3d does not support MI yet.

[Xinqiu Yao]

Q&A not a bug

[Rupesh Agarwal]

Program1:

library('bio3d')

library("MASS") trj <- read.dcd("production_1-50_every10.dcd")

pdb <- read.pdb("rxr_tr_ligs.pdb")

ca.inds <- atom.select(pdb, elety="CA")

xyz <- fit.xyz(fixed=pdb$xyz, mobile=trj, fixed.inds=ca.inds$xyz, mobile.inds=ca.inds$xyz)

cij<-dccm(xyz[,ca.inds$xyz], ncore=32, method="pearson")

write.matrix(cij, file="cross_corelation_dccm_pearson.dat", sep = " ")

and Program 2:

library('bio3d')

library("MASS")

trj <- read.dcd("production_1-50_every10.dcd")

pdb <- read.pdb("rxr_tr_ligs.pdb")

ca.inds <- atom.select(pdb, elety="CA")

xyz <- fit.xyz(fixed=pdb$xyz, mobile=trj, fixed.inds=ca.inds$xyz, mobile.inds=ca.inds$xyz)

cij<-dccm(xyz[,ca.inds$xyz], ncore=32, method="lmi")

write.matrix(cij, file="cross_corelation_dccm_lmi.dat", sep = " ")

These two programs give identical results. thoughts?

[Xinqiu Yao]

Could you tell us what version of bio3d you are using?

Note that the option method='lmi' etc. is not supported in the officially released version (on CRAN) yet. You need to install the development version. See here for more detail.