4 Copyright (C) 2003-2009 Kevin Thornton, krthornt[]@[]uci.edu
6 Remove the brackets to email me.
8 This file is part of libsequence.
10 libsequence is free software: you can redistribute it and/or modify
11 it under the terms of the GNU General Public License as published by
12 the Free Software Foundation, either version 3 of the License, or
13 (at your option) any later version.
15 libsequence is distributed in the hope that it will be useful,
16 but WITHOUT ANY WARRANTY; without even the implied warranty of
17 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18 GNU General Public License for more details.
20 You should have received a copy of the GNU General Public License
21 long with libsequence. If not, see <http://www.gnu.org/licenses/>.
25 #ifndef __SEQUENCE_BITS_SNN_TCC__
26 #define __SEQUENCE_BITS_SNN_TCC__
30 template< typename shuffler >
31 std::pair<double,double>
32 Snn_test_details(const PolyTable & snpTable,
33 const unsigned config[],
36 const unsigned & nperms)
38 std::vector < std::vector<double> > dkj(snpTable.size(),
39 std::vector<double>(snpTable.size(),0.));
40 for(unsigned i=0;i<snpTable.size()-1;++i)
42 for(unsigned j=i+1;j<snpTable.size();++j)
44 dkj[i][j] = double(NumDiffs(snpTable[i],snpTable[j]));
47 std::vector <unsigned> __individuals(snpTable.size(),0u);
48 for(unsigned i=0;i<snpTable.size();++i)
52 const double observed = Snn_statistic(&__individuals[0],
55 npop,unsigned(__individuals.size()));
56 unsigned perm = nperms;
60 //std::random_shuffle(__individuals.begin(),__individuals.end(),uni_int);
61 s(__individuals.begin(),__individuals.end());
62 double permuted_stat = Snn_statistic(&__individuals[0],
65 npop,unsigned(__individuals.size()));
66 pv += (permuted_stat >= observed) ? 1 : 0;
68 return std::make_pair(observed,double(pv)/double(nperms));
71 template< typename shuffler >
72 std::vector< std::vector<double> >
73 Snn_test_pairwise_details(const PolyTable & snpTable,
74 const unsigned config[],
77 const unsigned & nperms)
79 std::vector<unsigned> offsets(npop+1);
80 for(unsigned i=1;i<=npop;++i)
82 offsets[i] = offsets[i-1]+config[i-1];
84 std::vector < std::vector<double> > dkj(snpTable.size(),
85 std::vector<double>(snpTable.size(),0.));
86 for(unsigned i=0;i<snpTable.size()-1;++i)
88 for(unsigned j=i+1;j<snpTable.size();++j)
90 dkj[i][j] = Sequence::NumDiffs(snpTable[i],snpTable[j]);
93 std::vector <unsigned> __individuals(snpTable.size(),0u),individuals;
96 for(unsigned i=0;i<snpTable.size();++i)
100 std::vector< std::vector<double> > rv;
101 for(unsigned i=0;i< npop-1;++i)
103 for(unsigned j=i+1;j<npop;++j)
105 individuals.assign(__individuals.begin()+offsets[i],
106 __individuals.begin()+offsets[i+1]);
107 individuals.insert(individuals.end(),
108 __individuals.begin()+offsets[j],
109 __individuals.begin()+offsets[j+1]);
110 _config[0] = config[i];
111 _config[1] = config[j];
112 const double observed = Snn_statistic( &individuals[0],dkj,_config,2,
113 _config[0]+_config[1]);
114 unsigned p =0, _nperms = nperms;
117 //std::random_shuffle(__individuals.begin(),__individuals.end(),uni_int);
118 s(__individuals.begin(),__individuals.end());
119 double perm = Snn_statistic( &individuals[0],dkj,_config,2,
120 _config[0]+_config[1]);
121 if(perm>=observed) ++p;
124 std::vector<double> result(4);
125 result[0]=double(i+1);
126 result[1]=double(j+1);
127 result[2]=double(observed);
128 result[3]=double(p)/double(nperms);
129 rv.push_back(result);
135 template< typename shuffler >
136 std::pair<double,double>
137 Snn_test(const PolyTable & snpTable,
138 const unsigned config[],
141 const unsigned & nperms)
143 @brief Conducts a permutation-test of Hudson's Snn (sequence nearest-neighbor) statistic
144 \param snpTable The data on which we wish to perform the test.
145 \param config An array of the sample sizes in each deme.
146 \param npop The number of populations. For example, npop could equal config.size() if config were a vector
147 \param uni_int A random number generator whose operator() takes one argument, n, and returns a value uniformly-distributed on the half-open interval [0,N)
148 \param nperms The number of permutations to do for the test
149 \return A pair of doubles (std::pair<double,double>). the first member of the pair is
150 the observed value of the statistic, and the second member is the estimated p-value
151 \ingroup popgenanalysis
154 return Snn_test_details(snpTable,config,npop,s,nperms);
157 template< typename shuffler >
158 std::vector< std::vector<double> >
159 Snn_test_pairwise(const PolyTable & snpTable,
160 const unsigned config[],
163 const unsigned & nperms)
165 @brief Conducts a permutation-test of Hudson's Snn (sequence nearest-neighbor) statistic,
166 for all pairwise combinations of populations
167 \param snpTable The data on which we wish to perform the test.
168 \param config An array of the sample sizes in each deme.
169 \param npop The number of populations. For example, npop could equal config.size() if config were a vector
170 \param uni_int A random number generator whose operator() takes one argument, n,
171 and returns a value uniformly-distributed on the half-open interval [0,N)
172 \param nperms The number of permutations to do for the test
173 \return A vector of vector<double>. Each vector contains 4 elements, indexed 0 to 3. Elements 0 and 1
174 are the indexes of the i-th and j-th population. Element 2 is the observed Snn between the i-th and j-th
175 population, and element 3 is the p-value estimated by permutation
176 \ingroup popgenanalysis
179 return Snn_test_pairwise_details(snpTable,config,npop,s,nperms);
182 // template< typename shuffler >
183 // std::pair<double,double>
184 // Snn_test(const PolyTable & snpTable,
185 // const unsigned config[],
186 // const size_t & npop,
188 // const unsigned & nperms)
190 // @brief Conducts a permutation-test of Hudson's Snn (sequence nearest-neighbor) statistic
191 // \param snpTable The data on which we wish to perform the test.
192 // \param config An array of the sample sizes in each deme.
193 // \param npop The number of populations. For example, npop could equal config.size() if config were a vector
194 // \param uni_int A random number generator whose operator() takes one argument, n, and returns a value uniformly-distributed on the half-open interval [0,N)
195 // \param nperms The number of permutations to do for the test
196 // \return A pair of doubles (std::pair<double,double>). the first member of the pair is
197 // the observed value of the statistic, and the second member is the estimated p-value
198 // \ingroup popgenanalysis
201 // return Snn_test_details(snpTable,config,npop,s,nperms);
204 // template< typename shuffler >
205 // std::vector< std::vector<double> >
206 // Snn_test_pairwise(const PolyTable & snpTable,
207 // const unsigned config[],
208 // const size_t & npop,
210 // const unsigned & nperms)
212 // @brief Conducts a permutation-test of Hudson's Snn (sequence nearest-neighbor) statistic,
213 // for all pairwise combinations of populations
214 // \param snpTable The data on which we wish to perform the test.
215 // \param config An array of the sample sizes in each deme.
216 // \param npop The number of populations. For example, npop could equal config.size() if config were a vector
217 // \param uni_int A random number generator whose operator() takes one argument, n,
218 // and returns a value uniformly-distributed on the half-open interval [0,N)
219 // \param nperms The number of permutations to do for the test
220 // \return A vector of vector<double>. Each vector contains 4 elements, indexed 0 to 3. Elements 0 and 1
221 // are the indexes of the i-th and j-th population. Element 2 is the observed Snn between the i-th and j-th
222 // population, and element 3 is the p-value estimated by permutation
223 // \ingroup popgenanalysis
226 // return Snn_test_pairwise_details(snpTable,config,npop,s,nperms);