ARB
pvp.cxx
Go to the documentation of this file.
1 // ============================================================= //
2 // //
3 // File : pvp.cxx //
4 // Purpose : calculate positional variability (parsimony) //
5 // //
6 // Institute of Microbiology (Technical University Munich) //
7 // http://www.arb-home.de/ //
8 // //
9 // ============================================================= //
10 
11 #include "pvp.h"
12 #include "AP_pos_var.h"
13 
14 #include <TreeNode.h>
15 #include <arb_progress.h>
16 
17 GB_ERROR PVP_calculate(GBDATA *gb_main, const char *ali_name, const char *tree_name, const char *target_SAI_name) {
19  GB_transaction ta(gb_main);
20  long ali_len = GBT_get_alignment_len(gb_main, ali_name);
21 
22  if (ali_len <= 0) {
23  error = GB_await_error();
24  }
25  else {
26  const double phase1_fraction = 1.0/650;
27  arb_progress progress(WEIGHTED, "Calculating positional variability", phase1_fraction);
28  progress.subtitle("Loading Tree");
29 
30  // get tree
31  TreeNode *tree;
32  {
33  tree = GBT_read_tree(gb_main, tree_name, new SimpleRoot);
34  if (!tree) {
35  error = GB_await_error();
36  }
37  else {
38  GBT_link_tree(tree, gb_main, true, NULp, NULp);
39  }
40  }
41 
42  progress.inc_and_check_user_abort(error);
43 
44  if (!error) {
45  progress.subtitle("Counting mutations");
46 
47  GB_alignment_type at = GBT_get_alignment_type(gb_main, ali_name);
48  bool isNUC = at==GB_AT_DNA || at==GB_AT_RNA;
49 
50  AP_pos_var pvp(gb_main, ali_name, ali_len, isNUC, tree_name);
51 
52  error = pvp.delete_aliEntry_from_SAI(target_SAI_name);
53  if (!error) error = pvp.retrieve(tree);
54  if (!error) error = pvp.save_aliEntry_to_SAI(target_SAI_name);
55 
56  progress.inc_and_check_user_abort(error);
57  }
58 
59  destroy(tree);
60  }
61 
62  ta.close(error);
63 
64  return error;
65 }
66 
67 // --------------------------------------------------------------------------------
68 
69 #ifdef UNIT_TESTS
70 #ifndef TEST_UNIT_H
71 #include <test_unit.h>
72 #endif
73 
74 static arb_test::match_expectation saidata_equal(GBDATA *gb_main, const char *sainame, const char *aliname, const char *expected_data) {
75  using namespace arb_test;
76 
77  GB_transaction ta(gb_main);
78 
79  GBDATA *gb_sai = GBT_find_SAI(gb_main, sainame);
80  GBDATA *gb_saiali = GB_entry(gb_sai, aliname);
81  GBDATA *gb_saidata = GB_entry(gb_saiali, "data");
82 
83  const char *saidata = GB_read_char_pntr(gb_saidata);
84 
85  return that(saidata).is_equal_to(expected_data);
86 }
87 static arb_test::match_expectation saitype_equal(GBDATA *gb_main, const char *sainame, const char *aliname, const char *expected_type) {
88  using namespace arb_test;
89 
90  GB_transaction ta(gb_main);
91 
92  GBDATA *gb_sai = GBT_find_SAI(gb_main, sainame);
93  GBDATA *gb_saiali = GB_entry(gb_sai, aliname);
94  GBDATA *gb_saitype = GB_entry(gb_saiali, "_TYPE");
95 
96  const char *saitype = GB_read_char_pntr(gb_saitype);
97 
98  return that(saitype).is_equal_to(expected_type);
99 }
100 
101 #define TEST_EXPECT_SAIDATA_EQUAL(expected) TEST_EXPECTATION(saidata_equal(gb_main,SAI_name,aliname,expected))
102 #define TEST_EXPECT_SAITYPE_EQUAL(expected) TEST_EXPECTATION(saitype_equal(gb_main,SAI_name,aliname,expected))
103 
104 void TEST_pvp() {
105  GB_shell shell;
106  const char *SAI_name = "POS_VAR_BY_PARSIMONY";
107 
108  { // nuc data:
109  GBDATA *gb_main = GB_open("TEST_nuc.arb", "rw");
110  const char *aliname = "ali_16s";
111  TEST_EXPECT_NO_ERROR(PVP_calculate(gb_main, aliname, "tree_nuc", SAI_name));
112 
113  TEST_EXPECT_SAITYPE_EQUAL("PVP: Positional Variability by Parsimony: tree 'tree_nuc' ntaxa 12");
114  TEST_EXPECT_SAIDATA_EQUAL(".----222----7------32-7-7----774-----77-----747-77-7-47-----------73--433-.3134--1-41-4-4011222337--"
115  "....-632320124--11141-4-..24312434--3----7-77---7----774-4-4-34--44-77-732777424--7-----3247------44"
116  "33----------2--73063442221633737777326110244-3----73--2--4-773--4-7--3774777137-7774437274-3372--7--"
117  "--7------47-7---24--73--------2-4--7374---74-734------43-777-47-7---3---3----------7-3---3----------"
118  "-----7-----4-7----7--4-4------3--3---3----4-4---4-47-------77---4-----447-7--4-7-----377--7-7434373-"
119  "-77--5412243422747--7-441244321337--777--3434333-77---777---------7------------7-----------7---77--4"
120  "-7----477----4---7-----------4-774-4---7--742237-7-777-7---4---327447----4----43477423-----447-7--77"
121  "21234--7---7477-77777-727777---7-3--4---77------7------7---4--4--4-4---7---7------------47727--47777"
122  "4477-----7--4-7747-------47---7-7-7--------------------------74--7----------43-7--7477-774734-6----7"
123  "47---6-42-447-774-4--------77-----7-34-------7-------747------77-7---------4------------4-----------"
124  "7-----------------------4--4---7-.---7--------74437-------7-3424-732233427--7-72333127-6.--7777---6-"
125  "-424347-------7-------47--------------77----77---------------7------------------77437-77------4-743-"
126  "77-3-77---7----7-3377-------7-7-74--7-7---------7-7---4-------74---7----7------7-3447----7----------"
127  "--------2443----7777-33--2-4444-7--44437---37---7746----34447-47---7----7--4--4----------4--4--7----"
128  "44-----7---------7-4777-----3-77774-7------------7----47---------------------7--7-47---47-774777744-"
129  "--34-7-44227----3...7712.1----123-7742-4-777774-7774----4-3--7--------------------77--447--744--77-7"
130  "---------------3-...");
131 
132  // tests reported errors:
133  TEST_EXPECT_ERROR_CONTAINS(PVP_calculate(gb_main, "ali_61s", "tree_nuc", SAI_name), "alignment 'ali_61s' not found");
134  TEST_EXPECT_ERROR_CONTAINS(PVP_calculate(gb_main, "ali_16s", "tree_miss", SAI_name), "tree not found");
135 
136  GB_close(gb_main);
137  }
138 
139  { // protein data:
140  GBDATA *gb_main = GB_open("TEST_prot.arb", "rw");
141  const char *aliname = "ali_tuf_pro";
142  TEST_EXPECT_NO_ERROR(PVP_calculate(gb_main, aliname, "tree_prot", SAI_name));
143 
144  TEST_EXPECT_SAITYPE_EQUAL("PVP: Positional Variability by Parsimony: tree 'tree_prot' ntaxa 11");
145  TEST_EXPECT_SAIDATA_EQUAL("-----4-634363-4-3-4666----6--6-6-6666466666664666466666-43-34366633334446-4666-364-6--66666664-44646"
146  "666-4---4-66---6--66-4-66-6666444-6666466666--4--4--6466633346-62--3-3436343243266464433446433434662"
147  "2-4634-663466644346314-23364666646243343-6366-34632-0632-46646-666-43-6-6---46-46-3-3432-32634464334"
148  "443636666-66461612636----64664642421634-3-33444443632334-1-3-64-44664646634-433366666623624426246444"
149  "4436330334446-63434344622344666320664664--666643--66-446341446--6------6--466----------");
150 
151  GB_close(gb_main);
152  }
153 }
154 
155 #endif // UNIT_TESTS
156 
157 // --------------------------------------------------------------------------------
158 
159 
160 
const char * GB_ERROR
Definition: arb_core.h:25
GBDATA * GB_open(const char *path, const char *opent)
Definition: ad_load.cxx:1363
TreeNode * GBT_read_tree(GBDATA *gb_main, const char *tree_name, TreeRoot *troot)
Definition: adtree.cxx:791
long GBT_get_alignment_len(GBDATA *gb_main, const char *aliname)
Definition: adali.cxx:716
GB_ERROR GBT_link_tree(TreeNode *tree, GBDATA *gb_main, bool show_status, int *zombies, int *duplicates)
Definition: adtree.cxx:907
GBDATA * GBT_find_SAI(GBDATA *gb_main, const char *name)
Definition: aditem.cxx:177
GB_ERROR GB_await_error()
Definition: arb_msg.cxx:341
GB_ERROR save_aliEntry_to_SAI(const char *sai_name)
Definition: AP_pos_var.cxx:215
GB_ERROR delete_aliEntry_from_SAI(const char *sai_name)
Definition: AP_pos_var.cxx:202
static void error(const char *msg)
Definition: mkptypes.cxx:96
#define that(thing)
Definition: test_unit.h:1032
GB_alignment_type GBT_get_alignment_type(GBDATA *gb_main, const char *aliname)
Definition: adali.cxx:750
GB_ERROR PVP_calculate(GBDATA *gb_main, const char *ali_name, const char *tree_name, const char *target_SAI_name)
Definition: pvp.cxx:17
static SearchTree * tree[SEARCH_PATTERNS]
Definition: ED4_search.cxx:629
GB_alignment_type
Definition: arbdb_base.h:61
GB_ERROR retrieve(TreeNode *tree)
Definition: AP_pos_var.cxx:125
GB_ERROR close(GB_ERROR error)
Definition: arbdbpp.cxx:32
void subtitle(const char *stitle)
Definition: arb_progress.h:321
#define TEST_EXPECT_NO_ERROR(call)
Definition: test_unit.h:1107
#define NULp
Definition: cxxforward.h:114
#define TEST_EXPECT_ERROR_CONTAINS(call, part)
Definition: test_unit.h:1103
GB_transaction ta(gb_var)
void destroy(TreeNode *that)
Definition: TreeNode.h:560
GB_CSTR GB_read_char_pntr(GBDATA *gbd)
Definition: arbdb.cxx:904
GBDATA * gb_main
Definition: adname.cxx:33
GBDATA * GB_entry(GBDATA *father, const char *key)
Definition: adquery.cxx:334
void inc_and_check_user_abort(GB_ERROR &error)
Definition: arb_progress.h:332
void GB_close(GBDATA *gbd)
Definition: arbdb.cxx:655