3 from __future__
import print_function
5 from IMP
import ArgumentParser
7 __doc__ =
"Score each of a set of combinations."
13 colors[
"Rpt1"] = [0.78, 0.78, 0.73]
14 colors[
"Rpt2"] = [0.78, 0.66, 0.58]
15 colors[
"Rpt3"] = [0.77, 0.43, 0.5]
16 colors[
"Rpt4"] = [0.76, 0.29, 0.67]
17 colors[
"Rpt5"] = [0.51, 0.14, 0.75]
18 colors[
"Rpt6"] = [0.0, 0., 0.75]
19 colors[
"Rpn1"] = [0.34, 0.36, 0.27]
20 colors[
"Rpn2"] = [0.42, 0.43, 0.36]
21 colors[
"Rpn3"] = [0.49, 0.5, 0.44]
22 colors[
"Rpn5"] = [0.56, 0.57, 0.51]
23 colors[
"Rpn6"] = [0.64, 0.64, 0.59]
24 colors[
"Rpn7"] = [0.71, 0.71, 0.66]
25 colors[
"Rpn8"] = [0.78, 0.78, 0.74]
26 colors[
"Rpn9"] = [1, 0, 0]
27 colors[
"Rpn10"] = [0, 1, 0]
28 colors[
"Rpn11"] = [0, 0, 1]
29 colors[
"Rpn12"] = [0.5, 0.2, 0.4]
30 colors[
"a1"] = [0.78, 0.78, 0.73]
31 colors[
"a2"] = [0.78, 0.66, 0.58]
32 colors[
"a3"] = [0.77, 0.43, 0.5]
33 colors[
"a4"] = [0.76, 0.29, 0.67]
34 colors[
"a5"] = [0.51, 0.14, 0.75]
35 colors[
"a6"] = [0.0, 0., 0.75]
36 colors[
"a7"] = [0.34, 0.36, 0.27]
37 colors[
"a8"] = [0.42, 0.43, 0.36]
38 colors[
"a9"] = [0.49, 0.5, 0.44]
39 colors[
"a10"] = [0.56, 0.57, 0.51]
41 colors[
"a11"] = [0.78, 0.78, 0.73]
42 colors[
"a12"] = [0.78, 0.66, 0.58]
43 colors[
"a13"] = [0.77, 0.43, 0.5]
44 colors[
"a14"] = [0.76, 0.29, 0.67]
45 colors[
"a15"] = [0.51, 0.14, 0.75]
46 colors[
"a16"] = [0.0, 0., 0.75]
47 colors[
"a17"] = [0.34, 0.36, 0.27]
48 colors[
"a18"] = [0.42, 0.43, 0.36]
49 colors[
"a19"] = [0.49, 0.5, 0.44]
50 colors[
"a20"] = [0.56, 0.57, 0.51]
52 colors[
"a21"] = [0.78, 0.78, 0.73]
53 colors[
"a22"] = [0.78, 0.66, 0.58]
54 colors[
"a23"] = [0.77, 0.43, 0.5]
55 colors[
"a24"] = [0.76, 0.29, 0.67]
56 colors[
"a25"] = [0.51, 0.14, 0.75]
57 colors[
"a26"] = [0.0, 0., 0.75]
58 colors[
"a27"] = [0.34, 0.36, 0.27]
59 colors[
"a28"] = [0.42, 0.43, 0.36]
60 colors[
"a29"] = [0.49, 0.5, 0.44]
61 colors[
"a30"] = [0.56, 0.57, 0.51]
65 def decompose(dmap, mhs):
70 full_sampled_map.set_particles(all_ps)
71 full_sampled_map.resample()
72 full_sampled_map.calcRMS()
74 dmap.get_number_of_voxels(
76 ).dmean * full_sampled_map.get_header(
79 lower = dmap.get_number_of_voxels(
81 ).rms * full_sampled_map.get_header(
83 norm_factors = [upper, lower]
84 print(
"===============my norm factors:", upper, lower)
88 def score_each_protein(dmap, mhs, sd):
89 norm_factors = decompose(dmap, mhs)
91 mdl = mhs[0].get_model()
92 for i
in range(len(mhs)):
96 mh_dmap.set_particles(leaves)
100 sd.get_component_header(i).get_transformations_fn())
102 for fit
in fits[:15]:
113 scores.append(mh_scores)
114 print(
"=====mol", i, mh_scores)
119 usage =
"""%prog [options] <asmb> <asmb.proteomics> <asmb.mapping>
120 <alignment.params> <combinations> <score combinations [output]>
122 Score each of a set of combinations.
124 p = ArgumentParser(usage)
125 p.add_argument(
"-m",
"--max", dest=
"max", type=int, default=999999999,
126 help=
"maximum number of fits considered")
127 p.add_argument(
"assembly_file", help=
"assembly file name")
128 p.add_argument(
"proteomics_file", help=
"proteomics file name")
129 p.add_argument(
"mapping_file", help=
"mapping file name")
130 p.add_argument(
"param_file", help=
"parameter file name")
131 p.add_argument(
"combinations_file", help=
"combinations file name")
132 p.add_argument(
"scores_file", help=
"output scores file name")
133 return p.parse_args()
136 def run(asmb_fn, proteomics_fn, mapping_fn, params_fn, combs_fn,
137 scored_comb_output_fn, max_comb):
139 dmap = IMP.em.read_map(asmb.get_assembly_header().get_dens_fn())
140 dmap.get_header().set_resolution(
142 dmap.update_voxel_size(asmb.get_assembly_header().get_spacing())
143 dmap.set_origin(asmb.get_assembly_header().get_origin())
144 threshold = asmb.get_assembly_header().get_threshold()
147 colors = get_color_map()
148 names = list(colors.keys())
150 alignment_params = IMP.multifit.AlignmentParams(params_fn)
151 alignment_params.show()
154 print(
"=========", combs_fn)
162 prot_data, mapping_fn)
164 em_anchors = mapping_data.get_anchors()
170 mapping_data, asmb, alignment_params)
171 align.set_fast_scoring(
False)
173 mdl = align.get_model()
174 mhs = align.get_molecules()
175 align.add_states_and_filters()
176 rbs = align.get_rigid_bodies()
178 align.set_density_map(dmap, threshold)
180 for i, mh
in enumerate(mhs):
181 ensmb.add_component_and_fits(mh,
184 rgb = colors[mh.get_name()]
186 rgb = colors[names[i]]
189 for p in IMP.core.get_leaves(mh):
190 g= IMP.display.XYZRGeometry(p)
198 align.add_all_restraints()
200 rs = align.get_restraint_set().get_restraints()
201 print(
"Get number of restraints:", len(rs))
203 rr = IMP.RestraintSet.get_from(r)
204 for i
in range(rr.get_number_of_restraints()):
205 print(rr.get_restraint(i).get_name())
206 output = open(scored_comb_output_fn,
"w")
210 for i
in range(asmb.get_number_of_component_headers()):
211 c = asmb.get_component_header(i)
212 fn = c.get_reference_fn()
217 rr = IMP.RestraintSet.get_from(r)
218 for i
in range(rr.get_number_of_restraints()):
219 output.write(rr.get_restraint(i).get_name() +
"|")
224 print(
"Number of combinations:", len(combs[:max_comb]))
226 print(
"native score")
229 rr = IMP.RestraintSet.get_from(r)
230 for j
in range(rr.get_number_of_restraints()):
231 print(rr.get_restraint(j).get_name(), rr.evaluate(
False))
234 for i, comb
in enumerate(combs[:max_comb]):
235 print(
"Scoring combination:", comb)
236 ensmb.load_combination(comb)
239 rr = IMP.RestraintSet.get_from(r)
240 for j
in range(rr.get_number_of_restraints()):
241 print(rr.get_restraint(j).get_name())
242 rscore = rr.evaluate(
False)
244 num_violated = num_violated + 1
246 print(str(all_leaves[0]) +
" :: " + str(all_leaves[-1]))
247 score = sf.evaluate(
False)
249 msg =
"COMB" + str(i) +
"|"
251 rr = IMP.RestraintSet.get_from(r)
252 for j
in range(rr.get_number_of_restraints()):
253 current_name = rr.get_restraint(j).get_name()
254 if current_name != prev_name:
255 msg +=
' ' + current_name +
' '
256 prev_name = current_name
257 rscore = rr.get_restraint(j).
evaluate(
False)
258 msg += str(rscore) +
"|"
260 num_violated = num_violated + 1
264 num_violated) +
"||||" + str(
265 fitr.evaluate(
False)) +
"||:"
269 output.write(msg +
"\n")
271 ensmb.unload_combination(comb)
277 run(args.assembly_file, args.proteomics_file, args.mapping_file,
278 args.param_file, args.combinations_file, args.scores_file, args.max)
280 if __name__ ==
"__main__":
An ensemble of fitting solutions.
double get_coarse_cc_coefficient(const DensityMap *grid1, const DensityMap *grid2, double grid2_voxel_data_threshold, bool allow_padding=false, FloatPair norm_factors=FloatPair(0., 0.))
Calculates the cross correlation coefficient between two maps.
void write_pdb(const Selection &mhd, TextOutput out, unsigned int model=1)
def evaluate
Evaluate the score of the restraint.
Create a scoring function on a list of restraints.
SettingsData * read_settings(const char *filename)
GenericHierarchies get_leaves(Hierarchy mhd)
Get all the leaves of the bit of hierarchy.
void read_pdb(TextInput input, int model, Hierarchy h)
ProteinsAnchorsSamplingSpace read_protein_anchors_mapping(multifit::ProteomicsData *prots, const std::string &anchors_prot_map_fn, int max_paths=INT_MAX)
Align proteomics graph to EM density map.
Class for sampling a density map from particles.
double get_rmsd(const Selection &s0, const Selection &s1)
void transform(XYZ a, const algebra::Transformation3D &tr)
Apply a transformation to the particle.
Fitting atomic structures into a cryo-electron microscopy density map.
ProteomicsData * read_proteomics_data(const char *proteomics_fn)
Proteomics reader.
void set_log_level(LogLevel l)
Set the current global log level.
IntsList read_paths(const char *txt_filename, int max_paths=INT_MAX)
Read paths.
Calculate score based on fit to EM map.
FittingSolutionRecords read_fitting_solutions(const char *fitting_fn)
Fitting solutions reader.
double get_resolution(Model *m, ParticleIndex pi)
Estimate the resolution of the hierarchy as used by Representation.