Hi Bianca,
If the domains are rotated between the trajectories, the motion cannot be
captured by a single linear component. Think of projecting an arc on a
straight line. You loose the (deflection) part on the second component.
That looks like squashing. I've mentioned this in relation to MD
simulations in a paper about TRAIL and DR5.
Hope it helps,
Tsjerk
On Tue, Jan 14, 2014 at 3:48 AM, Bianca Ihrig <bian...@bii.a-star.edu.sg>wrote:
> Hi everyone,
>
> I have been doing PCA using bio3d. My protein of interest has two domains
> and there are two crystal structures, one in a closed (domains are close
> together) and one in an open state (domains are nearly opposite each
> other). I ran MD simulations on both states and did PCA. When doing PCA on
> either of the trajectories it is working fine but when I combine the two
> trajectories and all the frames to only one of the domains (as I would like
> to see the motion of the other domain for closed vs. open state).
> But when doing so the second domain (the one I am not align to) is
> "squashed" when viewing the vectors in pymol. I dont understand why and I
> dont know how to fix it.
> Many thanks in advance for your time and help!
>
> This is the command I am using:
>
> R
> library(bio3d)
>
> #read pdb file for open structure
> pdbo <- read.pdb("ref-open(fr1800).pdb")
> #to check: print(pdbo)
> #read pdb file for closed structure
> pdbc <- read.pdb("ref-closed(fr9990).pdb")
>
> #read traj files
> tc <- read.ncdf("A2-closed-1001fr-100ns.nc")
> to <- read.ncdf("A2-open-100ns-1000fr.nc")
>
> #combine the two traj (to and tc) into one called t2
> t2 <- rbind(tc[2:1001,],to[2:1001,]) # 2000 frames 10413 coordinates
>
> #Select residues of 1st domain (WHA) for alignment
> whacac <- atom.select(pdbc, resno="9,10,11,12,13,14,15,16,
> 17,18,19,20,21,22,23,24,25,26,27,47,48,49,50,51,52,53,54,56,
> 55,57,58,59,60,61,62,63,64,66,67,68,69,70,78,79,80,81,82", elety = "CA")
> #Select CA of the two pdb files
> cao <- atom.select(pdbo, elety = "CA") # 215 atoms
> cac <- atom.select(pdbc, elety = "CA") # 215 atoms
>
> # align the combined trajs to 1st domain of closed structure
> xyz <- fit.xyz(fixed = pdbc$xyz, mobile = t2, fixed.inds = whacac$xyz,
> mobile.inds = whacac$xyz)
> #to check use command: write.ncdf()
>
> #####-- Do PCA over CA of the whole protein (not only domain 1)
> pc <- pca.xyz(xyz[, cac$xyz])
> png("pca_open-closed.png")
> pc1 <- plot(pc)
> dev.off()
> #view pc vector in pymol
> view.modes(pc, mode = 1, outprefix = "pc1", scale = 30, dual = FALSE,
> launch = TRUE)
>
> Many thanks! Your help is much appreciated!
>
> ______________________________________________
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/
> posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
--
Tsjerk A. Wassenaar, Ph.D.
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
