haddock_restraints/core/commands/
ti.rs

1use crate::*;
2use std::error::Error;
3
4/// Analyzes the true interface of a protein structure and generates Ambiguous Interaction Restraints (AIRs).
5///
6/// This function reads a PDB file, identifies the true interface between chains based on a distance cutoff,
7/// creates interactors for each chain involved in the interface, and generates AIRs.
8///
9/// # Arguments
10///
11/// * `input_file` - A string slice that holds the path to the input PDB file.
12/// * `cutoff` - A reference to a f64 value specifying the distance cutoff (in Angstroms) for determining interfaces.
13///
14/// # Returns
15///
16/// A `Result<String, Box<dyn Error>>` which is Ok(String) if the function completes successfully, or an Error if something goes wrong.
17///
18pub fn true_interface(
19    pdb: pdbtbx::PDB,
20    cutoff: &f64,
21    pml: &Option<String>,
22) -> Result<String, Box<dyn Error>> {
23    let true_interface = get_true_interface(&pdb, *cutoff);
24    let chains_in_contact = get_chains_in_contact(&pdb, *cutoff);
25
26    // Sort the true_interface by chain id
27    let mut true_interface: Vec<_> = true_interface.iter().collect();
28    true_interface.sort_by(|a, b| a.0.cmp(b.0));
29
30    // NOTE: Here the IDs of the interactors are their position in the PDB file; this is so that
31    // we can handle any order of chains.
32    let mut interactors: Vec<Interactor> = Vec::new();
33    for (chain_id, residues) in true_interface.iter() {
34        // Get what is the position of this chain in the PDB, this will be its ID
35        let target_id = pdb
36            .chains()
37            .position(|chain| chain.id() == *chain_id)
38            .unwrap();
39
40        let mut interactor = Interactor::new(target_id as u16);
41        interactor.set_chain(chain_id);
42        interactor.set_active(residues.iter().map(|&residue| residue as i16).collect());
43        interactor.set_pdb(pdb.clone());
44        interactor.remove_buried_residues();
45
46        // Assign the targets
47        for (chain_i, chain_j) in chains_in_contact.iter() {
48            let target_chain = if chain_i == *chain_id {
49                chain_j
50            } else {
51                chain_i
52            };
53            if let Some(target_index) = pdb.chains().position(|chain| chain.id() == *target_chain) {
54                interactor.add_target(target_index as u16);
55            }
56        }
57
58        interactors.push(interactor);
59    }
60
61    // Make the restraints
62    let air = Air::new(interactors);
63    let tbl = air.gen_tbl().unwrap();
64
65    println!("{}", tbl);
66
67    if let Some(output_f) = pml {
68        air.gen_pml(output_f)
69    };
70
71    Ok(tbl)
72}
73
74/// Generates Unambiguous Topological Interactions (TIs) from a protein structure.
75///
76/// This function reads a PDB file, identifies the closest residue pairs based on a specified distance cutoff,
77/// and creates unambiguous interactors for each residue pair.
78///
79/// # Arguments
80///
81/// * input_file - A string slice that holds the path to the input PDB file.
82/// * cutoff - A reference to a f64 value specifying the distance cutoff (in Angstroms) for determining interactions.
83///
84/// # Returns
85///
86/// A Result<String, Box<dyn Error>> containing the generated TBL (Topological Restraints List) if successful.
87///
88pub fn unambig_ti(
89    pdb: pdbtbx::PDB,
90    cutoff: &f64,
91    pml: &Option<String>,
92) -> Result<String, Box<dyn Error>> {
93    let pairs = get_closest_residue_pairs(&pdb, *cutoff);
94
95    let mut interactors: Vec<Interactor> = Vec::new();
96    let mut counter = 0;
97    pairs.iter().for_each(|g| {
98        let mut interactor_i = Interactor::new(counter);
99        counter += 1;
100        let mut interactor_j = Interactor::new(counter);
101        interactor_j.add_target(counter - 1);
102        interactor_i.add_target(counter);
103        counter += 1;
104
105        interactor_i.set_chain(g.chain_i.as_str());
106        interactor_i.set_active(vec![g.res_i as i16]);
107        interactor_i.set_active_atoms(vec![g.atom_i.clone()]);
108        interactor_i.set_target_distance(g.distance);
109
110        interactor_j.set_chain(g.chain_j.as_str());
111        interactor_j.set_passive(vec![g.res_j as i16]);
112        interactor_j.set_passive_atoms(vec![g.atom_j.clone()]);
113
114        interactors.push(interactor_i);
115        interactors.push(interactor_j);
116    });
117
118    // Make the restraints
119    let air = Air::new(interactors);
120    let tbl = air.gen_tbl().unwrap();
121
122    println!("{}", tbl);
123
124    if let Some(output_f) = pml {
125        air.gen_pml(output_f)
126    };
127
128    Ok(tbl)
129}
130
131/// Lists the interface residues for each chain in a protein structure.
132///
133/// This function analyzes a PDB file to identify the interface residues between chains
134/// based on a specified distance cutoff, and prints the results.
135///
136/// # Arguments
137///
138/// * `input_file` - A string slice that holds the path to the input PDB file.
139/// * `cutoff` - A reference to a f64 value specifying the distance cutoff (in Angstroms)
140///   for determining interface residues.
141///
142/// # Returns
143///
144/// A `Result<(), Box<dyn Error>>` which is Ok(()) if the function completes successfully,
145/// or an Error if something goes wrong.
146///
147pub fn list_interface(pdb: pdbtbx::PDB, cutoff: &f64) -> Result<(), Box<dyn Error>> {
148    let true_interface = get_true_interface(&pdb, *cutoff);
149
150    for (chain_id, residues) in true_interface.iter() {
151        let mut sorted_res = residues.iter().collect::<Vec<_>>();
152        sorted_res.sort();
153
154        println!("Chain {}: {:?}", chain_id, sorted_res);
155    }
156
157    Ok(())
158}
159#[cfg(test)]
160mod tests {
161
162    use super::*;
163    use std::io::{BufReader, Cursor};
164
165    #[test]
166    fn test_true_interface() {
167        let expected_tbl = r#"assign ( resid 933 and segid A )
168       (
169        ( resid 46 and segid B )
170     or
171        ( resid 47 and segid B )
172       ) 2.0 2.0 0.0
173
174assign ( resid 950 and segid A )
175       (
176        ( resid 46 and segid B )
177     or
178        ( resid 47 and segid B )
179       ) 2.0 2.0 0.0
180
181assign ( resid 46 and segid B )
182       (
183        ( resid 933 and segid A )
184     or
185        ( resid 950 and segid A )
186       ) 2.0 2.0 0.0
187
188assign ( resid 47 and segid B )
189       (
190        ( resid 933 and segid A )
191     or
192        ( resid 950 and segid A )
193       ) 2.0 2.0 0.0
194
195"#;
196
197        let content = std::fs::read_to_string("tests/data/complex.pdb").unwrap();
198        let mut opts = pdbtbx::ReadOptions::new();
199        opts.set_format(pdbtbx::Format::Pdb)
200            .set_level(pdbtbx::StrictnessLevel::Loose);
201        let cursor = Cursor::new(content.into_bytes());
202        let reader = BufReader::new(cursor);
203        let (pdb, _) = opts.read_raw(reader).unwrap();
204
205        let opt: Option<String> = None;
206        match true_interface(pdb, &3.0, &opt) {
207            Ok(tbl) => assert_eq!(tbl, expected_tbl),
208            Err(_e) => (),
209        };
210    }
211    #[test]
212    fn test_true_interface_ba() {
213        let expected_tbl = r#"assign ( resid 933 and segid A )
214       (
215        ( resid 46 and segid B )
216     or
217        ( resid 47 and segid B )
218       ) 2.0 2.0 0.0
219
220assign ( resid 950 and segid A )
221       (
222        ( resid 46 and segid B )
223     or
224        ( resid 47 and segid B )
225       ) 2.0 2.0 0.0
226
227assign ( resid 46 and segid B )
228       (
229        ( resid 933 and segid A )
230     or
231        ( resid 950 and segid A )
232       ) 2.0 2.0 0.0
233
234assign ( resid 47 and segid B )
235       (
236        ( resid 933 and segid A )
237     or
238        ( resid 950 and segid A )
239       ) 2.0 2.0 0.0
240
241"#;
242
243        let content = std::fs::read_to_string("tests/data/complex_BA.pdb").unwrap();
244        let mut opts = pdbtbx::ReadOptions::new();
245        opts.set_format(pdbtbx::Format::Pdb)
246            .set_level(pdbtbx::StrictnessLevel::Loose);
247        let cursor = Cursor::new(content.into_bytes());
248        let reader = BufReader::new(cursor);
249        let (pdb, _) = opts.read_raw(reader).unwrap();
250
251        let opt: Option<String> = None;
252        match true_interface(pdb, &3.0, &opt) {
253            Ok(tbl) => assert_eq!(tbl, expected_tbl),
254            Err(_e) => (),
255        };
256    }
257
258    #[test]
259    fn test_unambigti() {
260        let content = std::fs::read_to_string("tests/data/two_res.pdb").unwrap();
261        let mut opts = pdbtbx::ReadOptions::new();
262        opts.set_format(pdbtbx::Format::Pdb)
263            .set_level(pdbtbx::StrictnessLevel::Loose);
264        let cursor = Cursor::new(content.into_bytes());
265        let reader = BufReader::new(cursor);
266        let (pdb, _) = opts.read_raw(reader).unwrap();
267
268        let opt: Option<String> = None;
269        let expected_tbl = "assign ( resid 2 and segid A and name CA ) ( resid 10 and segid B and name CA ) 9.1 2.0 0.0\n\n";
270
271        match unambig_ti(pdb, &5.0, &opt) {
272            Ok(tbl) => assert_eq!(tbl, expected_tbl),
273            Err(_e) => (),
274        }
275    }
276}
haddock_restraints/core/commands/ti.rs

haddock_restraints/core/commands/
ti.rs
