diff --git a/modules/base/pymod/CMakeLists.txt b/modules/base/pymod/CMakeLists.txt index 4eb36a9147d7d83b081d05f4abeb3ca0a2868467..22700288ae061f64b45fb741e13b5044b2bf8977 100644 --- a/modules/base/pymod/CMakeLists.txt +++ b/modules/base/pymod/CMakeLists.txt @@ -7,4 +7,4 @@ set(OST_BASE_PYMOD_SOURCES pymod(NAME base OUTPUT_DIR ost CPP ${OST_BASE_PYMOD_SOURCES} - PY __init__.py settings.py) + PY __init__.py settings.py stutil.py) diff --git a/modules/base/pymod/__init__.py b/modules/base/pymod/__init__.py index 3d2137f7e76f6a66a8cf2041319d801742fafedf..aa723917bec0aed0d90e3c682ff200c39e90b72a 100644 --- a/modules/base/pymod/__init__.py +++ b/modules/base/pymod/__init__.py @@ -17,6 +17,7 @@ # 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA #------------------------------------------------------------------------------ from _base import * +from stutil import * from ost import geom from ost import io diff --git a/modules/base/pymod/stutil.py b/modules/base/pymod/stutil.py new file mode 100644 index 0000000000000000000000000000000000000000..2f64beaefa3731c5689c73c585cbc3edb37d32f1 --- /dev/null +++ b/modules/base/pymod/stutil.py @@ -0,0 +1,95 @@ +import math +from ost import mol + +def FloatValueExtract(func): + """ + Decorator to wrap functions that take a list of float values. In addition to + passing in a list of float values directly, it is possible to extract the + values from attributes or generic properties. + """ + def _dec(xs, prop=None, attr=None): + if prop!=None: + if len(xs)==0: + return func([]) + assert attr==None + level=mol.Prop.Level.UNSPECIFIED + if isinstance(xs[0], mol.AtomBase): + level=mol.Prop.Level.ATOM + elif isinstance(xs[0], mol.ResidueBase): + level=mol.Prop.Level.RESIDUE + elif isinstance(xs[0], mol.ChainBase): + level=mol.Prop.Level.CHAIN + epm=mol.EntityPropertyMapper(prop, level) + vals=[] + for x in xs: + try: + vals.append(epm.Get(x)) + except: + pass + return func(vals) + if attr!=None: + vals=[] + for x in xs: + try: + vals.append(getattr(x, attr)) + except: + pass + return func(vals) + return func(xs) + return _dec + +@FloatValueExtract +def Mean(xs): + """ + Calculate mean of dataset + """ + if len(xs)==0: + raise RuntimeError("Can't calculate mean of empty sequence") + return sum(xs)/len(xs) + +@FloatValueExtract +def StdDev(xs): + """ + Calculate standard-deviation of dataset + + | sum[xi-<x>]^2 | + sigma=sqrt|---------------| + | n | + """ + mean=Mean(xs) + return math.sqrt(sum([(x-mean)**2 for x in xs])/len(xs)) + +@FloatValueExtract +def Min(xs): + return min(xs) + +@FloatValueExtract +def Max(xs): + return max(xs) + +def Correl(xs, ys): + """ + Calculates the correlation coefficient between xs and ys as + + sum[(xi-<x>)*(yi-<y>)] + r=---------------------- + (n-1)*sx*sy + + where <x>, <y> are the mean of dataset xs and ys, and, sx and sy are the + standard deviations. + """ + if len(xs)!=len(ys): + raise RuntimeError("Can't calculate correl. Sequence lengths do not match.") + if len(xs)==1: + raise RuntimeError("Can't calculate correl of sequences with length 1.") + mean_x=Mean(xs) + mean_y=Mean(ys) + sigma_x, sigma_y=(0.0, 0.0) + cross_term=0.0 + for x, y in zip(xs, ys): + cross_term+=(x-mean_x)*(y-mean_y) + sigma_x+=(x-mean_x)**2 + sigma_y+=(y-mean_y)**2 + sigma_x=math.sqrt(sigma_x/len(xs)) + sigma_y=math.sqrt(sigma_y/len(ys)) + return cross_term/((len(xs)-1)*sigma_x*sigma_y) diff --git a/modules/mol/base/pymod/CMakeLists.txt b/modules/mol/base/pymod/CMakeLists.txt index 590af690986b9c300a28917bdbc0aea48f6f0cd9..29c336dca643d9bbafac80ad9fb5e3959e6ef48c 100644 --- a/modules/mol/base/pymod/CMakeLists.txt +++ b/modules/mol/base/pymod/CMakeLists.txt @@ -18,6 +18,7 @@ export_query_view_wrapper.cc export_torsion.cc export_visitor.cc wrap_mol.cc +export_entity_property_mapper.cc ) pymod(NAME mol CPP ${OST_BASE_PYMOD_SOURCES} PY __init__.py) diff --git a/modules/mol/base/pymod/export_entity_property_mapper.cc b/modules/mol/base/pymod/export_entity_property_mapper.cc new file mode 100644 index 0000000000000000000000000000000000000000..fa271bed83a04df8773f8110e773b94e5a7676e1 --- /dev/null +++ b/modules/mol/base/pymod/export_entity_property_mapper.cc @@ -0,0 +1,68 @@ +#include <boost/python.hpp> + +#include <ost/mol/entity_property_mapper.hh> + +using namespace ost; +using namespace ost::mol; +using namespace boost::python; + +namespace { + + +Real (EntityPropertyMapper::*get_ah_a)(const AtomHandle&) const=&EntityPropertyMapper::Get; +Real (EntityPropertyMapper::*get_av_a)(const AtomView&) const=&EntityPropertyMapper::Get; +Real (EntityPropertyMapper::*get_rh_a)(const ResidueHandle&) const=&EntityPropertyMapper::Get; +Real (EntityPropertyMapper::*get_rv_a)(const ResidueView&) const=&EntityPropertyMapper::Get; +Real (EntityPropertyMapper::*get_ch_a)(const ChainHandle&) const=&EntityPropertyMapper::Get; +Real (EntityPropertyMapper::*get_cv_a)(const ChainView&) const=&EntityPropertyMapper::Get; + + +Real (EntityPropertyMapper::*get_ah_b)(const AtomHandle&, Real) const=&EntityPropertyMapper::Get; +Real (EntityPropertyMapper::*get_av_b)(const AtomView&, Real) const=&EntityPropertyMapper::Get; +Real (EntityPropertyMapper::*get_rh_b)(const ResidueHandle&, Real) const=&EntityPropertyMapper::Get; +Real (EntityPropertyMapper::*get_rv_b)(const ResidueView&, Real) const=&EntityPropertyMapper::Get; +Real (EntityPropertyMapper::*get_ch_b)(const ChainHandle&, Real) const=&EntityPropertyMapper::Get; +Real (EntityPropertyMapper::*get_cv_b)(const ChainView&, Real) const=&EntityPropertyMapper::Get; + +EntityPropertyMapper create_epm(const String& prop_name, char level) +{ + switch(level) { + case 'A': + case 'a': + return EntityPropertyMapper(prop_name, Prop::ATOM); + case 'R': + case 'r': + return EntityPropertyMapper(prop_name, Prop::RESIDUE); + case 'C': + case 'c': + return EntityPropertyMapper(prop_name, Prop::CHAIN); + case 'U': + case 'u': + return EntityPropertyMapper(prop_name, Prop::UNSPECIFIED); + default: + throw std::runtime_error(String("unknown property level '")+level+"'"); + } +} + +} + +void export_EntityPropertyMapper() +{ + class_<EntityPropertyMapper>("EntityPropertyMapper", + init<const String&, + Prop::Level>(arg("level")=Prop::UNSPECIFIED)) + .def("__init__", &create_epm) + .def("Get", get_ah_a) + .def("Get", get_av_a) + .def("Get", get_rh_a) + .def("Get", get_rv_a) + .def("Get", get_ch_a) + .def("Get", get_cv_a) + .def("Get", get_ah_b) + .def("Get", get_av_b) + .def("Get", get_rh_b) + .def("Get", get_rv_b) + .def("Get", get_ch_b) + .def("Get", get_cv_b) + ; +} diff --git a/modules/mol/base/pymod/wrap_mol.cc b/modules/mol/base/pymod/wrap_mol.cc index 82a449f681fbe98866c42acb0b8581fefd66d331..ce55b9942f8f5fd378ce11c78b924b6b5691a44d 100644 --- a/modules/mol/base/pymod/wrap_mol.cc +++ b/modules/mol/base/pymod/wrap_mol.cc @@ -41,6 +41,7 @@ void export_CoordGroup(); void export_PropertyID(); void export_BoundingBox(); void export_QueryViewWrapper(); +void export_EntityPropertyMapper(); BOOST_PYTHON_MODULE(_mol) { export_Entity(); @@ -61,7 +62,7 @@ BOOST_PYTHON_MODULE(_mol) export_PropertyID(); export_BoundingBox(); export_QueryViewWrapper(); - + export_EntityPropertyMapper(); class_<Transform>("Transform", init<>()) .def("GetMatrix",&Transform::GetMatrix) .def("GetTransposedMatrix",&Transform::GetTransposedMatrix)