1 from __future__
import print_function
6 from IMP
import ArgumentParser
13 Create a GMM from either density file (.mrc), a pdb file (.pdb)
14 Will detect input format from extension.
15 Outputs as text and optionally as a density map
18 p = ArgumentParser(description=desc)
20 p.add_argument(
"-t",
"--covar_type", dest=
"covar_type", default=
'full',
21 choices=[
'spherical',
'tied',
'diag',
'full'],
22 help=
"covariance type for the GMM")
23 p.add_argument(
"-m",
"--out_map", dest=
"out_map", default=
'',
24 help=
"write out the gmm to an mrc file")
25 p.add_argument(
"-a",
"--apix", dest=
"apix", default=1.0, type=float,
26 help=
"if you don't provide a map, set the voxel_size "
27 "here (for sampling)")
28 p.add_argument(
"-n",
"--num_samples", dest=
"num_samples", default=1000000,
29 type=int, help=
"num samples to draw from the density map")
30 p.add_argument(
"-i",
"--num_iter", dest=
"num_iter", default=100, type=int,
31 help=
"num iterations of GMM")
32 p.add_argument(
"-s",
"--threshold", dest=
"threshold", default=0.0,
33 type=float, help=
"threshold for the map before sampling")
35 p.add_argument(
"-f",
"--force_radii", dest=
"force_radii", default=-1.0,
37 help=
"force radii to be this value (spherical) "
38 "-1 means deactivated ")
39 p.add_argument(
"-w",
"--force_weight", dest=
"force_weight", default=-1.0,
41 help=
"force weight to be this value (spherical) "
42 "-1 means deactivated ")
43 p.add_argument(
"-e",
"--force_weight_frac", dest=
"force_weight_frac",
44 action=
"store_true", default=
False,
45 help=
"force weight to be 1.0/(num centers). "
46 "Takes precedence over -w")
47 p.add_argument(
"-d",
"--use_dirichlet", dest=
"use_dirichlet",
48 default=
False, action=
"store_true",
49 help=
"use dirichlet process for fit")
51 p.add_argument(
"-k",
"--multiply_by_mass", dest=
"multiply_by_mass",
52 default=
False, action=
"store_true",
53 help=
"if set, will multiply all weights by the total mass "
54 "of the particles (PDB ONLY)")
55 p.add_argument(
"-x",
"--chain", dest=
"chain", default=
None,
56 help=
"If you passed a PDB file, read this chain")
58 p.add_argument(
"-z",
"--use_cpp", dest=
"use_cpp", default=
False,
60 help=
"EXPERIMENTAL. Uses the IMP GMM code. "
62 p.add_argument(
"data_file", help=
"data file name")
63 p.add_argument(
"n_centers", type=int, help=
"number of centers")
64 p.add_argument(
"out_file", help=
"output file name")
69 data_fn = args.data_file
70 ncenters = args.n_centers
71 out_txt_fn = args.out_file
74 if not os.path.isfile(data_fn):
75 raise Exception(
"The data file you entered: " + data_fn
79 ext = data_fn.split(
'.')[-1]
86 mh, chain=args.chain).get_selected_particles()
90 if args.multiply_by_mass:
98 dmap.set_was_used(
True)
99 print(
'sampling points')
100 pts = IMP.isd.sample_points_from_density(
101 dmap, args.num_samples, args.threshold)
103 raise ValueError(
"data_fn extension must be pdb or mrc")
110 if args.force_weight_frac:
111 force_weight = 1.0/ncenters
113 force_weight = args.force_weight
114 if force_weight != -1:
115 print(
'weight forced to', force_weight)
116 if not args.use_dirichlet:
118 pts, ncenters, mdl, density_ps, args.num_iter, args.covar_type,
119 force_radii=args.force_radii, force_weight=args.force_weight,
120 mass_multiplier=mass_multiplier)
123 pts, ncenters, mdl, density_ps, args.num_iter, args.covar_type,
124 mass_multiplier=mass_multiplier)
130 print(
"This option is experimental, only works if you "
133 density_ps = IMP.isd_emxl.fit_gaussians_to_density(
134 mdl, dmap, args.num_samples, ncenters, args.num_iter,
135 args.threshold, gmm_threshold)
138 comments = [
'Created by create_gmm.py, IMP.isd version %s'
141 comments.append(
'ncenters: %d' % ncenters)
142 for key
in (
'covar_type',
'apix',
'num_samples',
'num_iter',
143 'threshold',
'force_radii',
'force_weight',
144 'force_weight_frac',
'use_dirichlet',
'multiply_by_mass',
146 comments.append(
'%s: %s' % (key, repr(getattr(args, key))))
148 if args.out_map !=
'':
158 if __name__ ==
"__main__":
Select non water and non hydrogen atoms.
double get_mass(ResidueType c)
Get the mass from the residue type.
GenericHierarchies get_leaves(Hierarchy mhd)
Get all the leaves of the bit of hierarchy.
void read_pdb(TextInput input, int model, Hierarchy h)
Class for storing model, its restraints, constraints, and particles.
std::string get_module_version()
Return the version of this module, as a string.
std::string get_relative_path(std::string base, std::string relative)
Return a path to a file relative to another file.
Basic utilities for handling cryo-electron microscopy 3D density maps.
A decorator for a particle with x,y,z coordinates.
algebra::BoundingBoxD< 3 > get_bounding_box(const DensityMap *m)
Select hierarchy particles identified by the biological name.
Inferential scoring building on methods developed as part of the Inferential Structure Determination ...