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