# -*- tab-width:2;indent-tabs-mode:t;show-trailing-whitespace:t;rm-trailing-spaces:t -*-
# vi: set ts=2 noet:
#
# (c) Copyright Rosetta Commons Member Institutions.
# (c) This file is part of the Rosetta software suite and is made available under license.
# (c) The Rosetta software is developed by the contributing members of the Rosetta Commons.
# (c) For more information, see http://www.rosettacommons.org. Questions about this can be
# (c) addressed to University of Washington UW TechTransfer, email: license@u.washington.edu.
library(ggplot2)
feature_analyses <- c(feature_analyses, methods::new("FeaturesAnalysis",
id = "hbond_scores_vs_score_bb_seqsep",
author = "Matthew O'Meara",
brief_description = "",
feature_reporter_dependencies = c("StructureFeatures", "ResidueFeatures", "HBondFeatures"),
run=function(self, sample_sources, output_dir, output_formats){
sele <- "
CREATE TEMPORARY TABLE ref_hbonds (
tag TEXT,
don_resNum INTEGER,
don_atmNum INTEGER,
acc_resNum INTEGER,
acc_atmNum INTEGER,
don_chem_type TEXT,
acc_chem_type TEXT,
energy REAL,
PRIMARY KEY(tag, don_resNum, don_atmNum, acc_resNum, acc_atmNum));
INSERT INTO ref_hbonds SELECT
struct.tag,
don.resNum AS don_resNum, don.atmNum AS don_atmNum,
acc.resNum AS acc_resNum, acc.atmNum AS acc_atmNum,
don.HBChemType AS don_chem_type,
acc.HBChemType AS acc_chem_type,
hb.energy AS energy
FROM
ref.structures AS struct,
ref.hbonds AS hb,
ref.hbond_sites AS don, hbond_sites AS acc,
ref.hbond_sites_pdb AS don_pdb, ref.hbond_sites_pdb AS acc_pdb
WHERE
hb.struct_id = struct.struct_id AND
don.struct_id = hb.struct_id AND don.site_id = hb.don_id AND
acc.struct_id = hb.struct_id AND acc.site_id = hb.acc_id AND
don_pdb.struct_id = hb.struct_id AND don_pdb.site_id = hb.don_id AND
don_pdb.heavy_atom_temperature < 30 AND
acc_pdb.struct_id = hb.struct_id AND acc_pdb.site_id = hb.acc_id AND
acc_pdb.heavy_atom_temperature < 30 AND
don.HBChemType == 'hbdon_PBA' AND acc.HBChemType == 'hbacc_PBA';
CREATE UNIQUE INDEX IF NOT EXISTS new.hbond_sites_struct_id_resNum_atmNum ON hbond_sites(
struct_id, resNum, atmNum);
CREATE UNIQUE INDEX IF NOT EXISTS new.hbonds_struct_id_acc_id_don_id ON hbonds(
struct_id, acc_id, don_id);
CREATE UNIQUE INDEX IF NOT EXISTS new.structures_tag ON structures( tag );
SELECT
ref_hb.*,
CASE ref_hb.acc_resNum - ref_hb.don_resNum
WHEN 2 THEN '2'
WHEN 3 THEN '3'
WHEN 4 THEN '4'
WHEN -2 THEN '-2'
WHEN -3 THEN '-3'
WHEN -4 THEN '-4'
ELSE 'long' END AS seq_sep,
new_hb.energy AS new_energy
FROM
ref_hbonds AS ref_hb,
new.structures AS new_struct,
new.hbonds AS new_hb,
new.hbond_sites AS new_don, new.hbond_sites AS new_acc
WHERE
new_struct.tag = ref_hb.tag AND
new_don.struct_id = new_struct.struct_id AND
new_don.resNum = ref_hb.don_resNum AND
new_don.atmNum = ref_hb.don_atmNum AND
new_acc.struct_id = new_struct.struct_id AND
new_acc.resNum = ref_hb.acc_resNum AND
new_acc.atmNum = ref_hb.acc_atmNum AND
new_hb.struct_id = new_struct.struct_id AND
new_hb.acc_id = new_acc.site_id AND
new_hb.don_id = new_don.site_id;"
f <- query_sample_sources_against_ref(sample_sources, sele)
sele <- "
DROP TABLE ref_hbonds;"
query_sample_sources(sample_sources, sele, warn_zero_rows=F)
f <- na.omit(f, method="r")
f$energy_weighted <- f$energy
f$new_energy_weighted <- f$new_energy
#hbond_sr_bb weighting
f <- transform(f,
energy_weighted = ifelse(
don_chem_type == 'hbdon_PBA' &
acc_chem_type == 'hbacc_PBA' &
acc_resNum - don_resNum < 5 &
acc_resNum - don_resNum > -5,
energy_weighted*.585, energy_weighted))
#hbond_lr_bb weighting
f <- transform(f,
energy_weighted = ifelse(
don_chem_type == 'hbdon_PBA' &
acc_chem_type == 'hbacc_PBA' &
acc_resNum - don_resNum >= 5 &
acc_resNum - don_resNum <= -5,
energy_weighted*1.17, energy_weighted ))
#hbond_bb_sc weighting
f <- transform(f,
energy_weighted = ifelse(
(don_chem_type == 'hbdon_PBA' & acc_chem_type != 'hbacc_PBA') |
(don_chem_type != 'hbdon_PBA' & acc_chem_type == 'hbacc_PBA'),
energy_weighted*1.17, energy_weighted))
#hbond_sc weighting
f <- transform(f,
energy_weighted = ifelse(
don_chem_type != 'hbdon_PBA' & acc_chem_type != 'hbacc_PBA',
energy_weighted*1.1, energy_weighted))
#hbond_sr_bb hbond_lr_bb hbond_bb_sc weighting
f <- transform(f,
new_energy_weighted = ifelse(
don_chem_type == 'hbdon_PBA' | acc_chem_type == 'hbacc_PBA',
new_energy_weighted*1.17, new_energy_weighted))
#hbond_sc weighting
f <- transform(f,
new_energy_weighted = ifelse(
don_chem_type != 'hbdon_PBA' & acc_chem_type != 'hbacc_PBA',
new_energy_weighted*1.1, new_energy_weighted))
diagonal <- data.frame(
x=seq(-1.755, 0, length.out=200),
y=seq(-1.755, 0, length.out=200))
plot_id <- "hbond_score_vs_score_bb_seqsep"
p <- ggplot(data=f, aes(x=energy_weighted, y=new_energy_weighted)) +
theme_bw() +
geom_point(size=.4) +
geom_line(data=diagonal, aes(x=x, y=y), color="darkgray") +
stat_density2d(size=.2) +
facet_wrap(~seq_sep) +
ggtitle(paste("Weighted HBond Energy Rosetta vs Sp2 Potential:") ) +
scale_x_continuous("Weighted Rosetta HBond Energy", limits=c(-1.755, 0), breaks=c(-1.5, -1, -.5)) +
scale_y_continuous("Weighted Sp2 Potential HBond Energy", limits=c(-3, 0), breaks=c(-3, -2, -1))
if(nrow(sample_sources) <= 3){
p <- p + theme(legend.position="bottom", legend.direction="horizontal")
}
save_plots(self, plot_id, sample_sources, output_dir, output_formats)
})) # end FeaturesAnalysis
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.