Q: How to calculate RMSF for specific residues

[Renhua Sun]

I am a beginner, even a beninner to use commend line.

Here I ask for help about How to calculate RMSF for specific residues, like Residues 90, 95 in Chain P and Residues 80-90 in Chain A.

Thank you very much!

Renhua

[Lars Skjærven]

Can you please elaborate? From what are you looking to calculate RMSF? Typically a trajectory from simulation or ensemble of protein structures.

[Renhua Sun]

Thank you for the reply. Yes, it is a trajectory file (*.dcd) from MD simulation. I could calculate RMSF for the whole object following http://thegrantlab.org/bio3d/tutorials/trajectory-analysis So, I raised the above question for specific residues :)

[Lars Skjærven]

Excellent. So now you have a vector of rmsf values and you just need to find the correct indices for the residues you're interested in (note: normally we'll just plot the whole thing and label by residue number). Having a PDB file with the same number of atoms (and order) as you trajectory is then needed.

> trtfile <- system.file("examples/hivp.dcd", package="bio3d")
> trj <- read.dcd(trtfile)

> pdbfile <- system.file("examples/hivp.pdb", package="bio3d")
> pdb <- read.pdb(pdbfile)

> rf <- rmsf(trj)

> sele <- atom.select(pdb, chain="A", resno=10)
> rf[ sele$atom ]
[1] 20.71261

Note that in this example both the PDB and DCD file contain only CA atoms. If you have all atoms you can select only these atoms you're interested in:

> sele <- atom.select(pdb, elety=c("CA", "C", "N", "O"))

# residue numbers to group by
> resno <- pdb$atom$resno[sele$atom]

# mean rmsf value of mainchain atoms of each residue
> rf <- rmsf(xyz[, sele$xyz], grpby=resno)

(see more here: http://thegrantlab.org/bio3d/html/rmsf.html)

You'll then have to find the correct indices of rf in which you are interested in.

[Renhua Sun]

Thank you very much!