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);
49 
50  bool isNUC = at==GB_AT_DNA || at==GB_AT_RNA;
51  AP_pos_var pvp(gb_main, ali_name, ali_len, isNUC, tree_name);
52 
53  error = pvp.delete_aliEntry_from_SAI(target_SAI_name);
54  if (!error) error = pvp.retrieve(tree);
55  if (!error) error = pvp.save_aliEntry_to_SAI(target_SAI_name);
56 
57  progress.inc_and_check_user_abort(error);
58  }
59 
60  destroy(tree);
61  }
62 
63  ta.close(error);
64 
65  return error;
66 }
67 
68 // --------------------------------------------------------------------------------
69 
70 #ifdef UNIT_TESTS
71 #ifndef TEST_UNIT_H
72 #include <test_unit.h>
73 #endif
74 
75 static arb_test::match_expectation saidata_equal(GBDATA *gb_main, const char *sainame, const char *aliname, const char *expected_data) {
76  using namespace arb_test;
77 
78  GB_transaction ta(gb_main);
79 
80  GBDATA *gb_sai = GBT_find_SAI(gb_main, sainame);
81  GBDATA *gb_saiali = GB_entry(gb_sai, aliname);
82  GBDATA *gb_saidata = GB_entry(gb_saiali, "data");
83 
84  const char *saidata = GB_read_char_pntr(gb_saidata);
85 
86  return that(saidata).is_equal_to(expected_data);
87 }
88 static arb_test::match_expectation saitype_equal(GBDATA *gb_main, const char *sainame, const char *aliname, const char *expected_type) {
89  using namespace arb_test;
90 
91  GB_transaction ta(gb_main);
92 
93  GBDATA *gb_sai = GBT_find_SAI(gb_main, sainame);
94  GBDATA *gb_saiali = GB_entry(gb_sai, aliname);
95  GBDATA *gb_saitype = GB_entry(gb_saiali, "_TYPE");
96 
97  const char *saitype = GB_read_char_pntr(gb_saitype);
98 
99  return that(saitype).is_equal_to(expected_type);
100 }
101 
102 #define TEST_EXPECT_SAIDATA_EQUAL(expected) TEST_EXPECTATION(saidata_equal(gb_main,SAI_name,aliname,expected))
103 #define TEST_EXPECT_SAITYPE_EQUAL(expected) TEST_EXPECTATION(saitype_equal(gb_main,SAI_name,aliname,expected))
104 
105 void TEST_pvp() {
106  GB_shell shell;
107  const char *SAI_name = "POS_VAR_BY_PARSIMONY";
108 
109  { // nuc data:
110  GBDATA *gb_main = GB_open("TEST_nuc.arb", "rw");
111  const char *aliname = "ali_16s";
112  TEST_EXPECT_NO_ERROR(PVP_calculate(gb_main, aliname, "tree_nuc", SAI_name));
113 
114  TEST_EXPECT_SAITYPE_EQUAL("PVP: Positional Variability by Parsimony: tree 'tree_nuc' ntaxa 12");
115  TEST_EXPECT_SAIDATA_EQUAL(".----222----7------32-7-7----774-----77-----747-77-7-47-----------73--433-.3134--1-41-4-4011222337--"
116  "....-632320124--11141-4-..24312434--3----7-77---7----774-4-4-34--44-77-732777424--7-----3247------44"
117  "33----------2--73063442221633737777326110244-3----73--2--4-773--4-7--3774777137-7774437274-3372--7--"
118  "--7------47-7---24--73--------2-4--7374---74-734------43-777-47-7---3---3----------7-3---3----------"
119  "-----7-----4-7----7--4-4------3--3---3----4-4---4-47-------77---4-----447-7--4-7-----377--7-7434373-"
120  "-77--5412243422747--7-441244321337--777--3434333-77---777---------7------------7-----------7---77--4"
121  "-7----477----4---7-----------4-774-4---7--742237-7-777-7---4---327447----4----43477423-----447-7--77"
122  "21234--7---7477-77777-727777---7-3--4---77------7------7---4--4--4-4---7---7------------47727--47777"
123  "4477-----7--4-7747-------47---7-7-7--------------------------74--7----------43-7--7477-774734-6----7"
124  "47---6-42-447-774-4--------77-----7-34-------7-------747------77-7---------4------------4-----------"
125  "7-----------------------4--4---7-.---7--------74437-------7-3424-732233427--7-72333127-6.--7777---6-"
126  "-424347-------7-------47--------------77----77---------------7------------------77437-77------4-743-"
127  "77-3-77---7----7-3377-------7-7-74--7-7---------7-7---4-------74---7----7------7-3447----7----------"
128  "--------2443----7777-33--2-4444-7--44437---37---7746----34447-47---7----7--4--4----------4--4--7----"
129  "44-----7---------7-4777-----3-77774-7------------7----47---------------------7--7-47---47-774777744-"
130  "--34-7-44227----3...7712.1----123-7742-4-777774-7774----4-3--7--------------------77--447--744--77-7"
131  "---------------3-...");
132 
133  // tests reported errors:
134  TEST_EXPECT_ERROR_CONTAINS(PVP_calculate(gb_main, "ali_61s", "tree_nuc", SAI_name), "alignment 'ali_61s' not found");
135  TEST_EXPECT_ERROR_CONTAINS(PVP_calculate(gb_main, "ali_16s", "tree_miss", SAI_name), "tree not found");
136 
137  GB_close(gb_main);
138  }
139 
140  { // protein data:
141  GBDATA *gb_main = GB_open("TEST_prot.arb", "rw");
142  const char *aliname = "ali_tuf_pro";
143  TEST_EXPECT_NO_ERROR(PVP_calculate(gb_main, aliname, "tree_prot", SAI_name));
144 
145  TEST_EXPECT_SAITYPE_EQUAL("PVP: Positional Variability by Parsimony: tree 'tree_prot' ntaxa 11");
146  TEST_EXPECT_SAIDATA_EQUAL("-----4-634363-4-3-4666----6--6-6-6666466666664666466666-43-34366633334446-4666-364-6--66666664-44646"
147  "666-4---4-66---6--66-4-66-6666444-6666466666--4--4--6466633346-62--3-3436343243266464433446433434662"
148  "2-4634-663466644346314-23364666646243343-6366-34632-0632-46646-666-43-6-6---46-46-3-3432-32634464334"
149  "443636666-66461612636----64664642421634-3-33444443632334-1-3-64-44664646634-433366666623624426246444"
150  "4436330334446-63434344622344666320664664--666643--66-446341446--6------6--466----------");
151 
152  GB_close(gb_main);
153  }
154 }
155 
156 #endif // UNIT_TESTS
157 
158 // --------------------------------------------------------------------------------
159 
160 
161 
#define arb_assert(cond)
Definition: arb_assert.h:245
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:837
long GBT_get_alignment_len(GBDATA *gb_main, const char *aliname)
Definition: adali.cxx:833
GB_ERROR GBT_link_tree(TreeNode *tree, GBDATA *gb_main, bool show_status, int *zombies, int *duplicates)
Definition: adtree.cxx:953
GBDATA * GBT_find_SAI(GBDATA *gb_main, const char *name)
Definition: aditem.cxx:177
GB_ERROR GB_await_error()
Definition: arb_msg.cxx:342
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:1043
GB_alignment_type GBT_get_alignment_type(GBDATA *gb_main, const char *aliname)
Definition: adali.cxx:873
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:35
void subtitle(const char *stitle)
Definition: arb_progress.h:321
#define TEST_EXPECT_NO_ERROR(call)
Definition: test_unit.h:1118
#define NULp
Definition: cxxforward.h:114
#define TEST_EXPECT_ERROR_CONTAINS(call, part)
Definition: test_unit.h:1114
GB_transaction ta(gb_var)
void destroy(TreeNode *that)
Definition: TreeNode.h:600
GB_CSTR GB_read_char_pntr(GBDATA *gbd)
Definition: arbdb.cxx:904
GBDATA * gb_main
Definition: adname.cxx:32
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