Back to Article
Dihedral angles
Download Source

Dihedral angles

Function to calculate chi angles for a single residue

In [1]:
res.chiR <- function(pdb.in, res){

  # trim pdb to residue
  x <- trim.pdb(pdb.in, resno = res)

  # number of frames
  frames <- nrow(x$xyz)

  # residue name
  residue <- unique(x$atom$resid)

  # define atom selection based on residue type
  if (residue %in% c('DG', 'DA')) {
    chi.atom <- atom.select(x, elety = c("O4'", "C1'", "N9", "C4"))
  } else {
    chi.atom <- atom.select(x, elety = c("O4'", "C1'", "N1", "C2"))
  }

  # trim pdb to chi atoms
  x <- trim.pdb(x, chi.atom)

  # check atom position so that dihedral angles are properly calculated
  if (residue %in% c('DG', 'DA')) {
    atom.match <- match(x$atom$elety, c("O4'", "C1'", "N9", "C4"))
  } else {
    atom.match <- match(x$atom$elety, c("O4'", "C1'", "N1", "C2"))
  }

  # calculate chi angles for all frames
  sapply(1:frames, function(i){

    xyz <-  x$xyz[i,]

    # if the atom order is not correct, reorder the atoms
    if (isFALSE(all(atom.match == 1:4, TRUE))) {
      xyz_reordered <- numeric() #initialized reordered coordinates
      for (i in atom.match) { #add the coordinates in the correct order
        start_index <- 3 * (i - 1) + 1
        xyz_reordered <- c(xyz_reordered, xyz[start_index:(start_index + 2)])
      }
      xyz <- xyz_reordered #replace the original coordinates with the reordered ones
    }

    # calculate chi angle
    xyz %>% torsion.xyz(.)

  }) %>%
    as.data.frame() %>%
    magrittr::set_colnames(c(paste0(residue, res)))
}

Function to calculate chi angles for a range of residues

In [2]:
chiR <- function(pdb.in, chain = NULL, res.start = NULL, res.end = NULL){

  # pull chain names if none supplied
  if (is.null(chain)) {
    chain <- unique(pdb.in$atom$chain)
  }

  # trim pdb to chain
  x <- trim.pdb(pdb.in, chain = chain)

  # pull residue min from selected chain if none supplied
  if (is.null(res.start)) {
    res.start <- min(x$atom$resno)
  }

  # pull residue max from selected chain if none supplied
  if (is.null(res.end)) {
    res.end <- max(x$atom$resno)
  }

  # calculate chi angles for all residues
  # using res.chiR function, which calculates chi angles for a single residue
  # across all frames
  sapply(res.start:res.end, function(i){
    res.chiR(x, i)
  }) %>%
    bind_rows() %>%
    mutate(frame = 1:nrow(.)) %>%
    pivot_longer(-frame, names_to = "residue", values_to = "chi")

}