3 from __future__
import print_function
5 from IMP
import OptionParser
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)
92 mdl = mhs[0].get_model()
93 for i
in range(len(mhs)):
97 mh_dmap.set_particles(leaves)
101 sd.get_component_header(i).get_transformations_fn())
103 for fit
in fits[:15]:
107 cc.cross_correlation_coefficient(
114 scores.append(mh_scores)
115 print(
"=====mol", i, mh_scores)
120 usage =
"""%prog [options] <asmb> <asmb.proteomics> <asmb.mapping>
121 <alignment.params> <combinatins> <score combinations [output]>
123 Score each of a set of combinations.
125 parser = OptionParser(usage)
126 parser.add_option(
"-m",
"--max", dest=
"max", default=999999999,
127 help=
"maximum number of fits considered")
128 (options, args) = parser.parse_args()
130 parser.error(
"incorrect number of arguments")
131 return [options, args]
134 def run(asmb_fn, proteomics_fn, mapping_fn, params_fn, combs_fn,
135 scored_comb_output_fn, max_comb):
137 dmap = IMP.em.read_map(asmb.get_assembly_header().get_dens_fn())
138 dmap.get_header().set_resolution(
140 dmap.update_voxel_size(asmb.get_assembly_header().get_spacing())
141 dmap.set_origin(asmb.get_assembly_header().get_origin())
142 threshold = asmb.get_assembly_header().get_threshold()
145 colors = get_color_map()
146 names = list(colors.keys())
148 alignment_params = IMP.multifit.AlignmentParams(params_fn)
149 alignment_params.show()
152 print(
"=========", combs_fn)
160 prot_data, mapping_fn)
162 em_anchors = mapping_data.get_anchors()
168 mapping_data, asmb, alignment_params)
169 align.set_fast_scoring(
False)
171 mdl = align.get_model()
172 mhs = align.get_molecules()
173 align.add_states_and_filters()
174 rbs = align.get_rigid_bodies()
176 align.set_density_map(dmap, threshold)
178 for i, mh
in enumerate(mhs):
179 ensmb.add_component_and_fits(mh,
182 rgb = colors[mh.get_name()]
184 rgb = colors[names[i]]
187 for p in IMP.core.get_leaves(mh):
188 g= IMP.display.XYZRGeometry(p)
196 align.add_all_restraints()
198 rs = align.get_restraint_set().get_restraints()
199 print(
"Get number of restraints:", len(rs))
201 rr = IMP.RestraintSet.get_from(r)
202 for i
in range(rr.get_number_of_restraints()):
203 print(rr.get_restraint(i).get_name())
204 output = open(scored_comb_output_fn,
"w")
208 for i
in range(asmb.get_number_of_component_headers()):
209 c = asmb.get_component_header(i)
210 fn = c.get_reference_fn()
215 rr = IMP.RestraintSet.get_from(r)
216 for i
in range(rr.get_number_of_restraints()):
217 output.write(rr.get_restraint(i).get_name() +
"|")
222 print(
"Number of combinations:", len(combs[:max_comb]))
224 print(
"native score")
227 rr = IMP.RestraintSet.get_from(r)
228 for j
in range(rr.get_number_of_restraints()):
229 print(rr.get_restraint(j).get_name(), rr.evaluate(
False))
232 for i, comb
in enumerate(combs[:max_comb]):
233 print(
"Scoring combination:", comb)
234 ensmb.load_combination(comb)
237 rr = IMP.RestraintSet.get_from(r)
238 for j
in range(rr.get_number_of_restraints()):
239 print(rr.get_restraint(j).get_name())
240 rscore = rr.evaluate(
False)
242 num_violated = num_violated + 1
244 print(str(all_leaves[0]) +
" :: " + str(all_leaves[-1]))
245 score = sf.evaluate(
False)
247 msg =
"COMB" + str(i) +
"|"
249 rr = IMP.RestraintSet.get_from(r)
250 for j
in range(rr.get_number_of_restraints()):
251 current_name = rr.get_restraint(j).get_name()
252 if current_name != prev_name:
253 msg +=
' ' + current_name +
' '
254 prev_name = current_name
255 rscore = rr.get_restraint(j).evaluate(
False)
256 msg += str(rscore) +
"|"
258 num_violated = num_violated + 1
262 num_violated) +
"||||" + str(
263 fitr.evaluate(
False)) +
"||:"
267 output.write(msg +
"\n")
269 ensmb.unload_combination(comb)
274 (options, args) = usage()
275 run(args[0], args[1], args[2], args[3], args[4], args[5], int(options.max))
277 if __name__ ==
"__main__":
An ensemble of fitting solutions.
void write_pdb(const Selection &mhd, TextOutput out, unsigned int model=1)
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)
Responsible for performing coarse fitting between two density objects.
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.