package HackaMol::X::Vina; #ABSTRACT: HackaMol extension for running Autodock Vina use Moose; use MooseX::StrictConstructor; use Moose::Util::TypeConstraints; use Math::Vector::Real; use MooseX::Types::Path::Tiny qw(AbsPath) ; use HackaMol; # for building molecules use File::chdir; use namespace::autoclean; use Carp; with qw(HackaMol::X::ExtensionRole); has $_ => ( is => 'rw', isa => AbsPath, predicate => "has_$_", required => 1, coerce => 1, ) foreach ( qw( receptor ligand ) ); has 'save_mol' => ( is => 'rw', isa => 'Bool', default => 0, ); has $_ => ( is => 'rw', isa => 'Num', predicate => "has_$_", ) foreach qw(center_x center_y center_z size_x size_y size_z); has 'num_modes' => ( is => 'rw', isa => 'Int', predicate => "has_num_modes", default => 1, lazy => 1, ); has $_ => ( is => 'rw', isa => 'Int', predicate => "has_$_", ) foreach qw(energy_range exhaustiveness seed cpu); has 'center' => ( is => 'rw', isa => 'Math::Vector::Real', predicate => "has_center", trigger => \&_set_center, ); has 'size' => ( is => 'rw', isa => 'Math::Vector::Real', predicate => "has_size", trigger => \&_set_size, ); sub BUILD { my $self = shift; if ( $self->has_scratch ) { $self->scratch->mkpath unless ( $self->scratch->exists ); } # build in some defaults $self->in_fn("conf.txt") unless ($self->has_in_fn); $self->exe("~/bin/vina") unless $self->has_exe; unless ( $self->has_out_fn ) { my $outlig = $self->ligand->basename; $outlig =~ s/\.pdbqt/\_out\.pdbqt/; $self->out_fn($outlig); } unless ( $self->has_command ) { my $cmd = $self->build_command; $self->command($cmd); } return; } sub _set_center { my ( $self, $center, $old_center ) = @_; $self->center_x( $center->[0] ); $self->center_y( $center->[1] ); $self->center_z( $center->[2] ); } sub _set_size { my ( $self, $size, $old_size ) = @_; $self->size_x( $size->[0] ); $self->size_y( $size->[1] ); $self->size_z( $size->[2] ); } #required methods sub build_command { my $self = shift; my $cmd; $cmd = $self->exe; $cmd .= " --config " . $self->in_fn->stringify; # we always capture output return $cmd; } sub _build_map_in { # this builds the default behavior, can be set anew via new return sub { return ( shift->write_input ) }; } sub _build_map_out { # this builds the default behavior, can be set anew via new my $sub_cr = sub { my $self = shift; my $qr = qr/^\s+\d+\s+(-*\d+\.\d)/; my ( $stdout, $sterr ) = $self->capture_sys_command; my @be = map { m/$qr/; $1 } grep { m/$qr/ } split( "\n", $stdout ); return (@be); }; return $sub_cr; } sub dock { my $self = shift; my $num_modes = shift; $self->num_modes($num_modes) if defined($num_modes); $self->map_input; return $self->map_output; } sub dock_mol { # want this to return configurations of the molecule my $self = shift; my $num_modes = shift; $self->num_modes($num_modes) if defined($num_modes); $self->map_input; local $CWD = $self->scratch if ( $self->has_scratch ); my @bes = $self->map_output; # this is fragile... broken if map_out changed... my $mol = HackaMol -> new(hush_read => 1) -> read_file_mol($self->out_fn->stringify); $mol->push_score(@bes); return ($mol); } sub write_input { my $self = shift; my $input; $input .= sprintf( "%-15s = %-55s\n", 'out', $self->out_fn->stringify ); $input .= sprintf( "%-15s = %-55s\n", 'log', $self->log_fn->stringify ) if $self->has_log_fn; foreach my $cond ( qw(receptor ligand cpu num_modes energy_range exhaustiveness seed)) { my $condition = "has_$cond"; $input .= sprintf( "%-15s = %-55s\n", $cond, $self->$cond ) if $self->$condition; } foreach my $metric (qw(center_x center_y center_z size_x size_y size_z)) { $input .= sprintf( "%-15s = %-55s\n", $metric, $self->$metric ); } $self->in_fn->spew($input); return ($input); } __PACKAGE__->meta->make_immutable; 1; __END__ =pod =head1 NAME HackaMol::X::Vina - HackaMol extension for running Autodock Vina =head1 VERSION version 0.00_5 =head1 SYNOPSIS use Modern::Perl; use HackaMol; use HackaMol::X::Vina; use Math::Vector::Real; my $receptor = "receptor.pdbqt"; my $ligand = "lig.pdbqt", my $rmol = HackaMol -> new( hush_read=>1 ) -> read_file_mol( $receptor ); my $lmol = HackaMol -> new( hush_read=>1 ) -> read_file_mol( $ligand ); my $fh = $lmol->print_pdb("lig_out.pdb"); my @centers = map {$_ -> xyz} grep {$_ -> name eq "OH" } grep {$_ -> resname eq "TYR"} $rmol -> all_atoms; foreach my $center ( @centers ){ my $vina = HackaMol::X::Vina -> new( receptor => $receptor, ligand => $ligand, center => $center, size => V( 20, 20, 20 ), cpu => 4, exhaustiveness => 12, exe => '~/bin/vina', scratch => 'tmp', ); my $mol = $vina->dock_mol(3); # fill mol with 3 binding configurations printf ("Score: %6.1f\n", $mol->get_score($_) ) foreach (0 .. $mol->tmax); $mol->print_pdb_ts([0 .. $mol->tmax], $fh); } $_->segid("hgca") foreach $rmol->all_atoms; #for vmd rendering cartoons.. etc $rmol->print_pdb("receptor.pdb"); =head1 DESCRIPTION HackaMol::X::Vina provides an interface to AutoDock Vina. This class does not include the AutoDock Vina program, which is L<released under a very permissive Apache license|http://vina.scripps.edu/manual.html#license>, with few restrictions on commercial or non-commercial use, or on the derivative works, such is this. Follow these L<instructions | http://vina.scripps.edu/manual.html#installation> to acquire the program. Most importantly, if you use this interface effectively, please be sure to cite AutoDock Vina in your work: O. Trott, A. J. Olson, AutoDock Vina: improving the speed and accuracy of docking with a new scoring function, efficient optimization and multithreading, Journal of Computational Chemistry 31 (2010) 455-461 Since HackaMol has no pdbqt writing capabilities (yet, HackaMol can read pdbqt files), the user is required to provide those files. This is still a work in progress and the API may still change. Documentation will improve as API gets more stable... comments welcome! The automated testing reported on metacpan will likely give a bunch of fails until I have time to figure out how to skip tests calling on the vina program to run. =head1 METHODS =head2 write_input This method takes no arguments; it returns, as a scalar, the input constructed from attributes. This method is called by map_input method via the map_in attribute to write the configuration file for running Vina. =head2 map_input provided by L<HackaMol::X::ExtensionRole>. Writes the configuration file for Vina. See dock and dock_mol methods. =head2 map_output provided by L<HackaMol::X::ExtensionRole>. By default, this method returns the docking scores as an array. =head2 dock_mol this method takes the number of binding modes (Integer) as an argument (Int). The argument is optional, and the num_modes attribute is rewritten if passed. This method calls the map_input and map_output methods for preparing and running Vina. It loads the resulting pdbqt and scores into a L<HackaMol::Molecule> object. The scores are stored into the score attribute provided by the L<HackaMol::QmMolRole>. See the synopsis for an example. =head2 dock this method is similar to dock_mol, but returns only the scores. =head1 ATTRIBUTES =head2 mol isa L<HackaMol::Molecule> object that is 'ro' and provided by L<HackaMol::X::ExtensionRole>. =head2 map_in map_out these attributes are 'ro' CodeRefs that can be adjusted in a given instance of a class. These are provided by L<HackaMol::X::ExtensionRole>. Setting the map_in and map_out attributes are for advanced use. Defaults are provided that are used in the map_input and map_output methods. =head2 receptor ligand these attributes are 'rw' and coerced into L<Path::Tiny> objects using the AbsPath type provided by L<MooseX::Types::Path::Tiny>. Thus, setting the receptor or ligand attributes with a string will store the entire path to the file, which is provided to Vina via the input configuration file. The receptor and ligand attributes typically point to pdbqt files used for running the docking calculations. =head2 save_mol this attribute isa 'Bool' that is 'rw'. =head2 center this attribute isa Math::Vector::Real object that is 'rw'. This attribute comes with a trigger that writes the center_x, center_y, and center_z attributes that are used in Vina configuration files. =head2 center_x center_y center_z this attribute isa Num that is 'rw'. These attributes provide the center for the box that (with size_x, size_y, size_z) define the docking space searched by Vina. Using the center attribute may be more convenient since it has the same type as the coordinates in atoms. See the synopsis. =head2 size_x size_y size_z this attribute isa Num that is 'rw'. These attributes provide the edgelengths of the the box that (with center_x, center_y, center_z) define the docking space searched by Vina. =head2 num_modes this attribute isa Int that is 'rw'. It provides the requested number of binding modes (ligand configurations) for Vina via the configuration file. Vina may return a smaller number of configurations depending on energy_range or other factors (that need documentation). =head2 energy_range this attribute isa Int that is 'rw'. In kcal/mol, provides a window for the number of configurations to return. =head2 exhaustiveness this attribute isa Int that is 'rw'. The higher the number the more time Vina will take looking for optimal docking configurations. =head2 cpu this attribute isa Int that is 'rw'. By default Vina will try to use all the cores available. Setting this attribute will limit the number of cores used by Vina. =head2 scratch this attribute isa L<Path::Tiny> that is 'ro' and provided by L<HackaMol::PathRole>. Setting this attribute return a Path::Tiny object with absolute path that will be created if needed and then used for all Vina calculations to be run. =head2 in_fn this attribute isa L<Path::Tiny> that is 'rw' and provided by L<HackaMol::PathRole>. The default is set to conf.txt when the object is built using the new method. If many instances of Vina will be running at the same time in the same directory, this conf.txt will need to be unique for each one!!! The same applies to out_fn which is described next. =head2 out_fn this attribute isa L<Path::Tiny> that is 'rw' and provided by L<HackaMol::PathRole>. The default is set to a value derived from the the basename of the ligand attribute. i.e. out_fn is set to lig_out.pdbqt from /some/big/path/lig.pdbqt. The Vina default behavior is to write to /some/big/path/lig_out.pdbqt, is usually not wanted (by me anyway); thus, the default is always set and written to the configuration file. If many instances of Vina will be running at the same time in the same directory, the output will need to be unique for each one as described above. =head1 EXTENDS =over 4 =item * L<Moose::Object> =back =head1 CONSUMES =over 4 =item * L<HackaMol::ExeRole> =item * L<HackaMol::ExeRole|HackaMol::PathRole> =item * L<HackaMol::PathRole> =item * L<HackaMol::X::ExtensionRole> =back =head1 AUTHOR Demian Riccardi <demianriccardi@gmail.com> =head1 COPYRIGHT AND LICENSE This software is copyright (c) 2014 by Demian Riccardi. This is free software; you can redistribute it and/or modify it under the same terms as the Perl 5 programming language system itself. =cut