~ubuntu-branches/ubuntu/oneiric/bioperl/oneiric

ok($hsp->query_string, 'LRVCGVANSKALLTNVHGLNLENWQEELAQAKEPF-NLGRLIRLVKEYHLLN----PVIVDCTSSQAVAD-QYADFLREGFHVVTPNKKANTSSMDYYHQLRYAAEKSRRKFLYDTNVGAGLPVIENLQNLLNAGDELMKFSGILSGSLSYIFGKLDE-GMSFSEATTLAREMGYTEPDPRDDLSGMDVARKLLILARET-GRELELADIEIEPVLPAEFNAEGDVAAFMANLSQLDDLFAARVAKARDEGKVLRYVGNIDEDGVCRVKIAEVDGNDPLFKVKNGENALAFYSHYYQPLPLVLRGYGAGNDVTAAGVFADLLRTLS');

100

110

ok($hsp->hit_string, 'GVVTGITDSREMLLSRIGLPLEIWKVALRDLEKPRKDLGKLDLTDDAFAVVDDPDIDVVVELTGGIEVARELYLDALEEGKHVVTANKALNASHGDEYLAL---AEKSGVDVLYEAAVAGGIPIIKTLRELLATGDRILKIEGIFNGTTNFILSEMDEKGLPFSDVLAEAQELGYTEADPRDDVEGIDAARKLAILARIAFGIELELDDVYVEGISPITAEDISSADEFGYTLKLLDEAMRQRVEDAESGGEVLRYPTLIPE-------------DHPLASVKGSDNAVAVEGEAYG--PLMFYGPGAGAEPTASAVVADIVRIAR');

101

111

ok($hsp->homology_string, ' V G+ +S+ +L + GL LE W+ L ++P +LG+L + +++ V+V+ T VA Y D L EG HVVT NK N S D Y L AEKS LY+ V G+P+I+ L+ LL GD ++K GI +G+ ++I ++DE G+ FS+ A+E+GYTE DPRDD+ G+D ARKL ILAR G ELEL D+ +E + P F L LD+ RV A G+VLRY I E + PL VK +NA+A Y PL+ G GAG + TA+ V AD++R ');

102

112

ok(join(' ', $hsp->seq_inds('query', 'gap',1)), '532 548-551 562 649 690');

103

113

# one more

104

114

$hit = $result->next_hit;

105

115

ok($hit);

106

116

107

117

while( $result = $searchio->next_result ) { ok($result); }

108

118

109

119

113

123

$result = $searchio->next_result;

114

124

115

125

ok($result->database_name, 'yeast.aa');

116

ok($result->query_name, 'gi|5763811|emb|CAB53164.1| putative transposase [Yersinia pestis]');

126

ok($result->query_name, 'gi|5763811|emb|CAB53164.1|');

127

ok($result->query_description, 'putative transposase [Yersinia pestis]');

128

ok($result->query_accession, 'CAB53164.1');

117

129

ok($result->query_length, 340);

118

130

119

131

$hit = $result->next_hit;

120

132

ok(! $hit);

121

133

134

$searchio = new Bio::SearchIO(-format => 'blastxml',

135

-file => Bio::Root::IO->catfile('t','data','mus.bls.xml'));

136

137

$result = $searchio->next_result;

138

139

ok($result->database_name,'Hs15_up1000');

140

ok($result->query_name,'NM_011441_up_1000_chr1_4505586_r');

141

ok($result->query_description,'chr1:4505586-4506585');

142

ok($result->query_accession,'NM_011441_up_1000_chr1_4505586_r');

143

ok($result->query_length,'1000');

144

$hit = $result->next_hit;

145

ok($hit->name,'NM_001938_up_1000_chr1_93161154_f');

146

ok($hit->description,'chr1:93161154-93162153');

147

ok($hit->accession,'3153');

148

ok($hit->length,'1000');

122

149

}

123

150

$searchio = new Bio::SearchIO ('-format' => 'blast',

124

151

'-file' => Bio::Root::IO->catfile('t','data','ecolitst.bls'));

130

157

ok($result->database_letters, 1358990);

131

158

132

159

ok($result->algorithm, 'BLASTP');

133

ok($result->algorithm_version, '2.1.3');

160

ok($result->algorithm_version, qr/^2\.1\.3/);

134

161

ok($result->query_name, qr/gi|1786183|gb|AAC73113.1| (AE000111) aspartokinase I,\s+homoserine dehydrogenase I [Escherichia coli]/);

135

162

ok($result->query_length, 820);

136

ok($result->get_statistic('kappa')== 0.041);

137

ok($result->get_statistic('lambda'), 0.267);

138

ok($result->get_statistic('entropy') == 0.14);

163

ok($result->get_statistic('kappa'), '0.135');

164

ok($result->get_statistic('kappa_gapped'), '0.0410');

165

ok($result->get_statistic('lambda'), '0.319');

166

ok($result->get_statistic('lambda_gapped'), '0.267');

167

ok($result->get_statistic('entropy'), '0.383');

168

ok($result->get_statistic('entropy_gapped'), '0.140');

169

139

170

ok($result->get_statistic('dbletters'), 1358990);

140

171

ok($result->get_statistic('dbentries'), 4289);

141

ok($result->get_statistic('hsplength'), 47);

172

ok($result->get_statistic('effective_hsplength'), 47);

142

173

ok($result->get_statistic('effectivespace'), 894675611);

143

174

ok($result->get_parameter('matrix'), 'BLOSUM62');

144

175

ok($result->get_parameter('gapopen'), 11);

145

176

ok($result->get_parameter('gapext'), 1);

146

147

my @valid = ( [ 'gb|AAC73113.1|', 820, 'AAC73113.1', '0.0', 1567],

148

[ 'gb|AAC76922.1|', 810, 'AAC76922.1', '1e-91', 332],

149

[ 'gb|AAC76994.1|', 449, 'AAC76994.1', '3e-47', 184]);

177

ok($result->get_statistic('S2'), '92');

178

ok($result->get_statistic('S2_bits'), '40.0');

179

ok($result->get_parameter('expect'), '1.0e-03');

180

ok($result->get_statistic('num_extensions'), '82424');

181

182

183

my @valid = ( [ 'gb|AAC73113.1|', 820, 'AAC73113', '0', 1567],

184

[ 'gb|AAC76922.1|', 810, 'AAC76922', '1e-91', 332],

185

[ 'gb|AAC76994.1|', 449, 'AAC76994', '3e-47', 184]);

150

186

my $count = 0;

151

187

while( $hit = $result->next_hit ) {

152

188

my $d = shift @valid;

154

190

ok($hit->name, shift @$d);

155

191

ok($hit->length, shift @$d);

156

192

ok($hit->accession, shift @$d);

157

ok($hit->significance, shift @$d );

193

ok(sprintf("%g",$hit->significance), sprintf("%g",shift @$d) );

158

194

ok($hit->raw_score, shift @$d );

159

195

160

196

if( $count == 0 ) {

164

200

ok($hsp->hit->start, 1);

165

201

ok($hsp->hit->end, 820);

166

202

ok($hsp->length('hsp'), 820);

167

203

ok($hsp->start('hit'), $hsp->hit->start);

204

ok($hsp->end('query'), $hsp->query->end);

205

ok($hsp->strand('sbjct'), $hsp->subject->strand);# alias for hit

168

206

ok($hsp->evalue == 0.0);

169

207

ok($hsp->score, 4058);

170

208

ok($hsp->bits,1567);

186

224

ok($result->database_letters, 1358990);

187

225

ok($result->database_entries, 4289);

188

226

ok($result->algorithm, 'BLASTP');

189

ok($result->algorithm_version, '2.0MP-WashU');

227

ok($result->algorithm_version, qr/^2\.0MP\-WashU/);

190

228

ok($result->query_name, qr/gi|1786183|gb|AAC73113.1| (AE000111) aspartokinase I,\s+homoserine dehydrogenase I [Escherichia coli]/);

191

229

ok($result->query_accession, 'AAC73113.1');

192

230

197

235

ok($result->get_statistic('dbletters'), 1358990);

198

236

ok($result->get_statistic('dbentries'), 4289);

199

237

ok($result->get_parameter('matrix'), 'BLOSUM62');

200

201

@valid = ( [ 'gb|AAC73113.1|', 820, 'AAC73113.1', '0.0', 4141],

202

[ 'gb|AAC76922.1|', 810, 'AAC76922.1', '3.1e-86', 844],

203

[ 'gb|AAC76994.1|', 449, 'AAC76994.1', '2.8e-47', 483]);

238

ok($result->get_statistic('Frame+0_lambda_used'), '0.319');

239

ok($result->get_statistic('Frame+0_kappa_used'), '0.136');

240

ok($result->get_statistic('Frame+0_entropy_used'), '0.384');

241

242

ok($result->get_statistic('Frame+0_lambda_computed'), '0.319');

243

ok($result->get_statistic('Frame+0_kappa_computed'), '0.136');

244

ok($result->get_statistic('Frame+0_entropy_computed'), '0.384');

245

246

ok($result->get_statistic('Frame+0_lambda_gapped'), '0.244');

247

ok($result->get_statistic('Frame+0_kappa_gapped'), '0.0300');

248

ok($result->get_statistic('Frame+0_entropy_gapped'), '0.180');

249

250

@valid = ( [ 'gb|AAC73113.1|', 820, 'AAC73113', '0', 4141],

251

[ 'gb|AAC76922.1|', 810, 'AAC76922', '3.1e-86', 844],

252

[ 'gb|AAC76994.1|', 449, 'AAC76994', '2.8e-47', 483]);

204

253

$count = 0;

205

254

while( $hit = $result->next_hit ) {

206

255

my $d = shift @valid;

208

257

ok($hit->name, shift @$d);

209

258

ok($hit->length, shift @$d);

210

259

ok($hit->accession, shift @$d);

211

ok($hit->significance, shift @$d );

260

ok(sprintf("%g",$hit->significance), sprintf("%g",shift @$d) );

212

261

ok($hit->raw_score, shift @$d );

213

262

214

263

if( $count == 0 ) {

241

290

ok($result->database_letters, 4662239);

242

291

ok($result->database_entries, 400);

243

292

ok($result->algorithm, 'TBLASTX');

244

ok($result->algorithm_version, '2.1.2');

245

ok($result->query_name, qr/HUMBETGLOA Human haplotype C4 beta-globin gene, complete cds./);

293

ok($result->algorithm_version, qr/^2\.1\.2/);

294

ok($result->query_name, 'HUMBETGLOA');

295

ok($result->query_description, 'Human haplotype C4 beta-globin gene, complete cds.');

246

296

ok($result->query_length, 3002);

247

297

ok($result->get_statistic('kappa'), 0.135);

248

298

ok($result->get_statistic('lambda'), 0.318);

250

300

ok($result->get_statistic('dbletters'), 4662239);

251

301

ok($result->get_statistic('dbentries'), 400);

252

302

ok($result->get_statistic('T'), 13);

303

ok($result->get_statistic('X1'), 16);

304

ok($result->get_statistic('X1_bits'), 7.3);

305

ok($result->get_statistic('X2'), 0);

306

ok($result->get_statistic('X2_bits'), '0.0');

307

ok($result->get_statistic('S1'), 41);

308

ok($result->get_statistic('S1_bits'), 21.7);

309

ok($result->get_statistic('S2'), 53);

310

ok($result->get_statistic('S2_bits'), 27.2);

311

253

312

ok($result->get_statistic('decayconst'), 0.1);

254

313

255

314

ok($result->get_parameter('matrix'), 'BLOSUM62');

272

331

ok($hsp->query->start, 1057);

273

332

ok($hsp->query->end, 1134);

274

333

ok($hsp->query->strand, 1);

334

ok($hsp->strand('query'), $hsp->query->strand);

275

335

ok($hsp->hit->end, 5893);

276

336

ok($hsp->hit->start, 5816);

277

337

ok($hsp->hit->strand, -1);

338

ok($hsp->strand('sbjct'), $hsp->subject->strand);

278

339

ok($hsp->length('hsp'), 26);

279

340

280

341

ok($hsp->evalue == 0.13);

289

350

ok($hsp->query_string, 'SAYWSIFPPLGCWWSTLGPRGSLSPL');

290

351

ok($hsp->hit_string, 'AAVWALFPPVGSQWGCLASQWRTSPL');

291

352

ok($hsp->homology_string, '+A W++FPP+G W L + SPL');

353

ok(join(' ', $hsp->seq_inds('query', 'nomatch',1)), '355 364 365 367 368 370 371 373-375');

292

354

}

293

355

}

294

356

last if( $count++ > @valid );

302

364

ok($result->database_entries, 657);

303

365

ok($result->algorithm, 'FASTN');

304

366

ok($result->algorithm_version, '3.3t08');

305

ok($result->query_name, qr/HUMBETGLOA Human haplotype C4 beta-globin gene, complete cds./);

367

ok($result->query_name, "HUMBETGLOA");

368

ok($result->query_description, "Human haplotype C4 beta-globin gene, complete cds.");

306

369

ok($result->query_length, 3002);

307

370

ok($result->get_parameter('gapopen'), -16);

308

371

ok($result->get_parameter('gapext'), -4);

312

375

ok($result->get_statistic('dbletters'), 112936249);

313

376

ok($result->get_statistic('dbentries'), 657);

314

377

315

@valid = ( [ 'BACR21I23', 73982, 'BACR21I23', '0.017', 44],

316

[ 'BACR40P19', 73982, 'BACR40P19', '0.017', 44],

317

[ 'BACR30L17', 32481, 'BACR30L17', '0.018', 44]);

378

@valid = ( [ 'BACR21I23', 73982, 'BACR21I23', '0.017', 44.2],

379

[ 'BACR40P19', 73982, 'BACR40P19', '0.017', 44.2],

380

[ 'BACR30L17', 32481, 'BACR30L17', '0.018', 44.1]);

318

381

$count = 0;

319

382

320

383

while( my $hit = $result->next_hit ) {

347

410

ok($hsp->gaps('hit'),1);

348

411

ok($hsp->query_string, 'GATTAAAACCTTCTGGTAAGAAAAGAAAAAATATATATATATATATATGTGTATATGTACACACATACATATACATATATATGCATTCATTTGTTGTTGTTTTTCTTAATTTGCTCATGCATGCTA----ATAAATTATGTCTAAAAATAGAAT---AAATACAAATCAATGTGCTCTGTGCATTA-GTTACTTATTAGGTTTTGGGAAACAAGAGGTAAAAAACTAGAGACCTCTTAATGCAGTCAAAAATACAAATAAATAAAAAGTCACTTACAACCCAAAGTGTGACTATCAATGGGGTAATCAGTGGTGTCAAATAGGAGGT');

349

412

ok($hsp->hit_string, 'GATGTCCTTGGTGGATTATGGTGTTAGGGTATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATAATATAATACAAAATATAATACAAAATATAATACAAAATATAATACAAAATATAATACAAAATATAATACAAAATATAATACAAAATATAATATAAAATATAATATAAAATATAATATAAAATAAAATATAAAATAAAATATAAAATAAAATATAAAATAAAATATAAAATAAAATAT-AATATAAAATATAAAATAAAATATAATATAAAATATAATATAAAATATAATATAAAATATAATATAAAATA');

350

ok($hsp->homology_string, ' :::::::::::::::::: : ::::: :: : : ::: ::::: :::::::: :: :: : : : : : : : :: : :: :: :: : ::: ::: :::::: ::: ::::: :: ::: : : : :: ::: : :: : : : : :: : :: : : :: : : :: : : ::: ::: :: ::::: ::: : : :: :: ::: : : : ::: :: ');

413

414

ok(join(' ', $hsp->seq_inds('query', 'nomatch',1)), '33 37 39 41 43 47-49 52 55 56 58 60 64 70 71 74 78 82 84 86 87 90-96 98 100 103 105 107 110-112 114 117 119 121-123 125 127-129 132 134 135 139-141 143 145-148 150-153 155 156 160 161 164 170 173 180-184 188 192 194 196-198 201 204 206-209 212 213 215 217 219 221 223-225 227 229 232 233 236 237 246 252 256 258 260 263 269 271');

415

ok(join(' ', $hsp->seq_inds('query', 'conserved',1)), '31 32 34-36 38 40 42 44-46 50 51 53 54 57 59 61-63 65-69 72 73 75-77 79-81 83 85 88 89 97 99 101 102 104 106 108 109 113 115 116 118 120 124 126 130 131 133 136-138 141 142 144 149 154 157-159 162 163 165-172 174-179 185-187 189-191 193-195 199 200 202 203 205 210 211 214 216 218 220 222 226 228 230 231 234 235 238-245 247-251 253-255 257 259 261 262 264-268 270 272-289');

416

# note: the reason this is not the same percent id above

417

# is we are calculating average percent id

418

ok(sprintf("%.2f",$hsp->get_aln->percentage_identity()), '59.30');

351

419

}

352

420

}

353

421

last if( $count++ > @valid );

372

440

ok($result->get_statistic('dbentries'), 4289);

373

441

374

442

375

@valid = ( [ 'gi|1787478|gb|AAC74309.1|', 512, 'AAC74309', 1.2, 29],

376

[ 'gi|1790635|gb|AAC77148.1|', 251, 'AAC77148', 2.1, 27],

377

[ 'gi|1786590|gb|AAC73494.1|', 94, 'AAC73494', 2.1, 26]);

443

@valid = ( [ 'gi|1787478|gb|AAC74309.1|', 512, 'AAC74309', 1.2, 29.2],

444

[ 'gi|1790635|gb|AAC77148.1|', 251, 'AAC77148', 2.1, 27.4],

445

[ 'gi|1786590|gb|AAC73494.1|', 94, 'AAC73494', 2.1, 25.9]);

378

446

$count = 0;

379

447

380

448

while( my $hit = $result->next_hit ) {

381

449

my $d = shift @valid;

382

383

450

ok($hit->name, shift @$d);

384

451

ok($hit->length, shift @$d);

385

452

ok($hit->accession, shift @$d);

406

473

ok($hsp->gaps('query'), 7);

407

474

ok($hsp->gaps, 49);

408

475

ok($hsp->query_string, 'NKEAIFTDDLPVADYLDDEFINSIPTAFDWRTRGAVTPVKNQGQCGSCWSFSTT-GNV----EGQHFISQNKLVSLSEQNLVDCDHECME-YEGEEACDEGCNGGLQPNAYNYIIKNGGIQTESSYPYTAETGTQCNFNSANIGAKISNFTMIPKNETVMAGYIVSTGP-LAIAADAVEWQFYIGGVFDIPCNPNSLDHGILIVGYSAKNTIFRKNMPYWIVKNSWGADWGEQGYIYLRRGKNTCGVSNFVSTSII');

409

ok($hsp->hit_string, 'MKIRSQVGMVLNLDKCIGCHTCSVTCKNVWTSREGVEYAWFNNVETKPGQGF-PTDWENQEKYKGGWI--RKINGKLQPRMGNRAMLLGKIFANPHLPGIDDYYEPFDFDYQNLHTAPEG----SKSQPIARPRSLITGERMAKIEKGPNWEDDLGGEFDKLAKDKNFDN-IQKAMYSQFENTFMMYLPRLCEHCLNPACVATCPSGAIYKREEDGIVLIDQDKCRGWRMCITGCPYKKIYFNWKSGKSEKCIFCYPRIEAGQPTVCSETC');

410

ok($hsp->homology_string, ' . :. : : : .: .: . :.: :: :: .. :.. . :.. : : .: :.: . :: ::: : . : : .. : . .:. :. . . :.. . . :: .: . .:. .:: :: . ...:. : . :: .. : .: ');

411

}

412

}

413

last if( $count++ > @valid );

414

}

476

ok($hsp->hit_string, (' 'x29).'MKIRSQVGMVLNLDKCIGCHTCSVTCKNVWTSREGVEYAWFNNVETKPGQGF-PTDWENQEKYKGGWI--RKINGKLQPRMGNRAMLLGKIFANPHLPGIDDYYEPFDFDYQNLHTAPEG----SKSQPIARPRSLITGERMAKIEKGPNWEDDLGGEFDKLAKDKNFDN-IQKAMYSQFENTFMMYLPRLCEHCLNPACVATCPSGAIYKREEDGIVLIDQDKCRGWRMCITGCPYKKIYFNWKSGKSEKCIFCYPRIEAGQPTVCSETC');

477

ok($hsp->homology_string, ' . :. : : : .: .: . :.: :: :: .. :.. . :.. : : .: :.: . :: ::: : . : : .. : . .:. :. . . :.. . . :: .: . .:. .:: :: . ...:. : . :: .. : .: '.' 'x60);

478

# note: the reason this is not the same percent id above

479

# is we are calculating average percent id

480

ok(sprintf("%.2f",$hsp->get_aln->percentage_identity()), 26.01);

481

}

482

}

483

last if( $count++ > @valid );

484

}

485

ok($result->hits, 8);

486

$searchio = new Bio::SearchIO(-format => 'fasta',

487

-file => 't/data/cysprot_vs_gadfly.FASTA');

488

$result = $searchio->next_result;

489

ok($result->database_name, qr/gadflypep2/);

490

ok($result->database_letters, 7177762);

491

ok($result->database_entries, 14334);

492

ok($result->algorithm, 'FASTP');

493

ok($result->algorithm_version, '3.3t08');

494

ok($result->query_name, 'cysprot.fa');

495

ok($result->query_length, 2385);

496

ok($result->get_parameter('gapopen'), -12);

497

ok($result->get_parameter('gapext'), -2);

498

ok($result->get_parameter('ktup'), 2);

499

ok($result->get_parameter('matrix'), 'BL50');

500

501

ok($result->get_statistic('lambda'), 0.1397);

502

ok($result->get_statistic('dbletters'), 7177762 );

503

ok($result->get_statistic('dbentries'), 14334);

504

505

506

@valid = ( [ 'Cp1|FBgn0013770|pp-CT20780|FBan0006692', 341,

507

'FBan0006692', '3.1e-59', 227.8],

508

[ 'CG11459|FBgn0037396|pp-CT28891|FBan0011459', 336,

509

'FBan0011459', '6.4e-41', 166.9],

510

[ 'CG4847|FBgn0034229|pp-CT15577|FBan0004847', 390,

511

'FBan0004847', '2.5e-40', 165.2]);

512

$count = 0;

513

514

while( my $hit = $result->next_hit ) {

515

my $d = shift @valid;

516

517

ok($hit->name, shift @$d);

518

ok($hit->length, shift @$d);

519

ok($hit->accession, shift @$d);

520

ok(sprintf("%g",$hit->significance), sprintf("%g",shift @$d) );

521

ok($hit->raw_score, shift @$d );

522

523

if( $count == 0 ) {

524

while( my $hsp = $hit->next_hsp ) {

525

ok($hsp->query->start, 1373);

526

ok($hsp->query->end, 1706);

527

ok($hsp->query->strand, 0);

528

ok($hsp->hit->start, 5);

529

ok($hsp->hit->end, 341);

530

ok($hsp->hit->strand, 0);

531

ok($hsp->length('hsp'), 345);

532

ok(sprintf("%g",$hsp->evalue), sprintf("%g",'3.1e-59') );

533

ok($hsp->score, 1170.6);

534

ok($hsp->bits,227.8);

535

ok(sprintf("%.2f",$hsp->percent_identity), 53.04);

536

ok(sprintf("%.4f",$hsp->frac_identical('query')), 0.5479);

537

ok(sprintf("%.4f",$hsp->frac_identical('hit')), '0.5430');

538

ok($hsp->query->frame(), 0);

539

ok($hsp->hit->frame(), 0);

540

ok($hsp->gaps('query'), 11);

541

ok($hsp->gaps, 194);

542

ok($hsp->hit_string, (' 'x26).'MRTAVLLPLLAL----LAVAQA-VSFADVVMEEWHTFKLEHRKNYQDETEERFRLKIFNENKHKIAKHNQRFAEGKVSFKLAVNKYADLLHHEFRQLMNGFNYTLHKQLRAADESFKGVTFISPAHVTLPKSVDWRTKGAVTAVKDQGHCGSCWAFSSTGALEGQHFRKSGVLVSLSEQNLVDCSTKYGNNGCNGGLMDNAFRYIKDNGGIDTEKSYPYEAIDDSCHFNKGTVGATDRGFTDIPQGDEKKMAEAVATVGPVSVAIDASHESFQFYSEGVYNEPQCDAQNLDHGVLVVGFGTDESGED---YWLVKNSWGTTWGDKGFIKMLRNKENQCGIASASSYPLV');

543

ok($hsp->query_string, 'SNWGNNGYFLIERGKNMCGLAACASYPIPQVMNPTLILAAFCLGIASATLTFDHSLEAQWTKWKAMHNRLY-GMNEEGWRRAVWEKNMKMIELHNQEYREGKHSFTMAMNAFGDMTSEEFRQVMNGFQ---NRKPR------KGKVFQEPLFYEAPRSVDWREKGYVTPVKNQGQCGSCWAFSATGALEGQMFRKTGRLISLSEQNLVDCSGPQGNEGCNGGLMDYAFQYVQDNGGLDSEESYPYEATEESCKYNPKYSVANDTGFVDIPK-QEKALMKAVATVGPISVAIDAGHESFLFYKEGIYFEPDCSSEDMDHGVLVVGYGFESTESDNNKYWLVKNSWGEEWGMGGYVKMAKDRRNHCGIASAASYPTVMTPLLLLAVLCLGTALATPKFDQTFNAQWHQWKSTHRRLYGTNEE');

544

# note: the reason this is not the same percent id above

545

# is we are calculating average percent id

546

ok(sprintf("%.2f",$hsp->get_aln->percentage_identity()), 56.13);

547

}

548

}

549

last if( $count++ > @valid );

550

}

551

ok($result->hits, 21);

415

552

416

553

# test on TFASTXY

417

554

$searchio = new Bio::SearchIO(-format => 'fasta',

434

571

ok($result->get_statistic('dbentries'), 9190);

435

572

436

573

437

@valid = ( [ 'NR_SC:SW-YNN2_YEAST', 1056, 'NR_SC:SW-YNN2_YEAST','1.6e-154', 547],

438

[ 'NR_SC:SW-MPCP_YEAST', 311, 'NR_SC:SW-MPCP_YEAST', '1.3e-25', 117],

439

[ 'NR_SC:SW-YEO3_YEAST', 300, 'NR_SC:SW-YEO3_YEAST', '5.7e-05', 48]);

574

@valid = ( [ 'NR_SC:SW-YNN2_YEAST', 1056, 'NR_SC:SW-YNN2_YEAST','1.6e-154', '547.0'],

575

[ 'NR_SC:SW-MPCP_YEAST', 311, 'NR_SC:SW-MPCP_YEAST', '1.3e-25', 117.1],

576

[ 'NR_SC:SW-YEO3_YEAST', 300, 'NR_SC:SW-YEO3_YEAST', '5.7e-05', 48.5]);

440

577

$count = 0;

441

578

442

579

while( my $hit = $result->next_hit ) {

445

582

ok($hit->name, shift @$d);

446

583

ok($hit->length, shift @$d);

447

584

ok($hit->accession, shift @$d);

448

ok($hit->significance, shift @$d );

585

ok(sprintf("%g",$hit->significance), sprintf("%g",shift @$d) );

449

586

ok($hit->raw_score, shift @$d );

450

587

451

588

if( $count == 0 ) {

458

595

ok($hsp->hit->strand, 0);

459

596

ok($hsp->length('hsp'), 1165);

460

597

461

ok($hsp->evalue == 1.6e-154);

598

ok(sprintf("%g",$hsp->evalue), sprintf("%g",'1.6e-154'));

462

599

ok($hsp->score, 2877.6);

463

600

ok($hsp->bits,'547.0');

464

601

ok(sprintf("%.2f",$hsp->percent_identity), 51.67);

470

607

ok($hsp->query_string, 'RKQLDPRIPALINNGVKANHRSFFVMVGDKGRDQVCPGMQAAMRFD*HRCR/LVNLHFLLSQARVSSRPSVLWCYKKD-LGFTT*VAASENLQQTIYFRPIATSHRKKREAKIKRDVKRGIRDANEQDPFELFVTVTDIRYTYYKDSAKILGQTFGMLVLQDYEAITPNLLARTIETVEGGGIVVLLLKTMSSLKQLYAMAM/DKL*CRDGVE*SDFS*LLI*DVHSRYRTDAHQFVQPRFNERFILSLGSNPDCLVLDDELNVLPLSKGKDIQIGKAGEEDDRGRKRKAEELKEMKENLEGVDIVGSLAKLAKTVDQAKAILTFVEAISEKNLSSTVALTAGRGRGKSAALGLAIGAALAHDYSNIFVTSPDPENLKTLFEFVFKALDALGYEEHIDYDVVQSTNPDFKKAIVRVNIFRGHRQTIQYISPEDSHVLGQAELVIIDEAAAIPLPLVRKLIGPYLVFMASTINGYEGTGRSLSIKLIQQLREQTRPSITKDSENAAASSAGSSSKAAAAGRSGAGLVRSLREIKLDEPIRYSPGDNVEKWLNNLLCLDATIVSK---SIQGCPHPSKCELYYVNRDTLFSYHPASEVFLQRMMALYVASHYKNSPNDLQMLSDAPAHHLFVLLPPIDEND-NTLPDPLVVLQVALEGNISREAILKEMAQSGMRSSGDMIPWIISTQFQDNDFATLSGARVVRIATHPDYARMGYGSRAMEALESFYNGTSYNFDDVPVDMGESFAD\VPRSDL*VTSFIPFPQNRTSTECVSQNANLQNDTIAIRDPSRMPPLLQRLSERKPETLDYLGVSFGLTRDLLRFWKKGGFTPLYASQKENALTGEYTFVMLKVLASAGGGGEWLGAFAQGMSCLLLQDEVHMGND*RL*TDFRQRFMNLLSYEAFKKFDASIALSILESTVPRNSPSPAP----KLLTNTELSSLLTPFDIKRLESYADSMLDYHVVLDLVPTIASLFFGKRLETS--LPPAQQAILLALGLQRKNVEALENELGITSTQTLALFGKVLRKMTKSLEDIRKASIASELP-----AEPTLAGRSANGSNKFVALQQTIEQDLADSAVQLNGEDDDASKKEQRELLNTLNMEEFAI-DQGGDWTEAEKQVERLASGKGGTRLSSTVSVKVDKLDD\AKRRRRRARMRVPRMRRR');

471

608

ok($hsp->hit_string, 'KKAIDSRIPSLIRNGVQTKQRSIFVIVGDRARNQ------------------LPNLHYLMMSADLKMNKSVLWAYKKKLLGFT--------------------SHRKKRENKIKKEIKRGTREVNEMDPFESFISNQNIRYVYYKESEKILGNTYGMCILQDFEALTPNLLARTIETVEGGGIVVILLKSMSSLKQLYTMTM-D--------------------VHARYRTEAHGDVVARFNERFILSLGSNPNCLVVDDELNVLPLSGAKNVKPLPPKEDDELPPKQL--ELQELKESLEDVQPAGSLVSLSKTVNQAHAILSFIDAISEKTLNFTVALTAGRGRGKSAALGISIAAAVSHGYSNIFVTSPSPENLKTLFEFIFKGFDALGYQEHIDYDIIQSTNPDFNKAIVRVDIKRDHRQTIQYIVPQDHQVLGQAELVVIDEAAAIPLPIVKNLLGPYLVFMASTINGYEGTGRSLSLKLIQQLRNQNNTSGRESTQTAVVSRDNKEKDSHLHSQS-----RQLREISLDEPIRYAPGDPIEKWLNKLLCLDVTLIKNPRFATRGTPHPSQCNLFVVNRDTLFSYHPVSENFLEKMMALYVSSHYKNSPNDLQLMSDAPAHKLFVLLPPIDPKDGGRIPDPLCVIQIALEGEISKESVRNSLSR-GQRAGGDLIPWLISQQFQDEEFASLSGARIVRIATNPEYASMGYGSRAIELLRDYFEGKF-------TDMSE---D-VRPKDYSI--------KRVSDKELAKT-NLLKDDVKLRDAKTLPPLLLKLSEQPPHYLHYLGVSYGLTQSLHKFWKNNSFVPVYLRQTANDLTGEHTCVMLNVLE--GRESNWLVEFAK---------------------DFRKRFLSLLSYD-FHKFTAVQALSVIESSKKAQDLSDDEKHDNKELTRTHLDDIFSPFDLKRLDSYSNNLLDYHVIGDMIPMLALLYFGDKMGDSVKLSSVQSAILLAIGLQRKNIDTIAKELNLPSNQTIAMFAKIMRKMSQYFRQLLSQSIEETLPNIKDDAIAEMDGEEIKNYNAAEALDQ-MEEDLEEAG----SEAVQAMREKQKELINSLNLDKYAINDNSEEWAESQKSLEIAAKAKGVVSLKTGKKRTTEKAED-IYRQEMKA-MKKPRKSKK');

472

609

ok($hsp->homology_string, '.: .: :::.:: :::....::.::.:::..:.: : :::.:. .: .. :::: ::: :::: ::::::: :::...::: :..::.:::: :.. .:::.:::.: ::::.:.:: .:::.::.:::::::::::::::::::.:::.::::::::.:.: : ::.::::.:: : ::::::::::::::.:::.:::::::::: .:... :.:. :. ::.:.::.:: :. .:::..:.:::.::.:::.:..:::::.:. :::::::::::::::::..:.::..: :::::::::.::::::::::.::..:::::.::::::..:::::::.::::::.: : :::::::: :.: .::::::::.::::::::::.:..:.::::::::::::::::::::::.:::::::.:. : .....:..: .. . . ..: :.::::.:::::::.::: .:::::.:::::.:.... . .: ::::.:.:. :::::::::::.:: ::..::::::.:::::::::::..::::::.::::::::: .: . .:::: :.:.::::.::.:.. . ... :.:..::.:::.:: ::::..::.:::::.:::::.:.:: :::::::.: :.....: .::.: : : .: . .:.: . .... :: .: . .:: . .:::: .:::. :. : :::::.:::..: .:::...:.:.: : : ::::.: :::.:: : ..:: ::. :::.::..::::. :.:: : :::..::. .. : : :: :.:.....:::.:::.::....:::::. :..: .: :.:: .. : : .:.:::::.::::::.... .::.. :.::.:.:.:..:::.. .... . :: :: : . :. .. : ::.: .:.:: ... .: .: ...:.::.:.::....:: :.. .:.:..:..: :..:: . :.. . ..: .: :.. .: :. :: ..');

610

# note: the reason this is not the same percent id above

611

# is we are calculating average percent id

612

ok(sprintf("%.2f",$hsp->get_aln->overall_percentage_identity()),

613

'51.77');

614

ok(sprintf("%.2f",$hsp->get_aln->average_percentage_identity()),

615

'58.41');

473

616

}

474

617

}

475

618

last if( $count++ > @valid );

476

619

}

477

620

ok($result->hits, 58);

478

621

# test for MarkW bug in blastN

479

622

480

623

$searchio = new Bio::SearchIO('-format' => 'blast',

486

629

ok($result->database_letters, 4677375331);

487

630

ok($result->database_entries, 1083200);

488

631

ok($result->algorithm, 'BLASTN');

489

ok($result->algorithm_version, '2.2.1');

632

ok($result->algorithm_version, qr/^2\.2\.1/);

490

633

ok($result->query_name, '');

491

634

ok($result->query_length, 60);

492

635

ok($result->get_parameter('gapopen'), 5);

498

641

ok($result->get_statistic('entropy'),1.31 );

499

642

ok($result->get_statistic('T'), 0);

500

643

ok($result->get_statistic('A'), 30);

501

ok($result->get_statistic('X1'), 6);

644

ok($result->get_statistic('X1'), '6');

645

ok($result->get_statistic('X1_bits'), 11.9);

502

646

ok($result->get_statistic('X2'), 15);

647

ok($result->get_statistic('X2_bits'), 29.7);

503

648

ok($result->get_statistic('S1'), 12);

649

ok($result->get_statistic('S1_bits'), 24.3);

504

650

ok($result->get_statistic('S2'), 17);

651

ok($result->get_statistic('S2_bits'), 34.2);

505

652

506

653

ok($result->get_statistic('dbentries'), 1083200);

507

654

508

@valid = ( [ 'gb|AY052359.1|', 2826, 'AY052359.1', '3e-18', 96],

509

[ 'gb|AC002329.2|AC002329', 76170, 'AC002329', '3e-18', 96],

510

[ 'gb|AF132318.1|AF132318', 5383, 'AF132318', '0.040', 42]);

655

@valid = ( [ 'gb|AY052359.1|', 2826, 'AY052359', '3e-18', 96, 1, 60,

656

'1.0000'],

657

[ 'gb|AC002329.2|AC002329', 76170, 'AC002329', '3e-18', 96, 1, 60,

658

'1.0000' ],

659

[ 'gb|AF132318.1|AF132318', 5383, 'AF132318', '0.04', 42, 35, 55,

660

'0.3500']);

511

661

$count = 0;

512

662

513

663

while( my $hit = $result->next_hit ) {

515

665

ok($hit->name, shift @$d);

516

666

ok($hit->length, shift @$d);

517

667

ok($hit->accession, shift @$d);

518

ok($hit->significance, shift @$d );

668

ok(sprintf("%g",$hit->significance), sprintf("%g",shift @$d) );

519

669

ok($hit->raw_score, shift @$d );

520

670

ok($hit->start, shift @$d);

671

ok($hit->end,shift @$d);

672

ok(sprintf("%.4f",$hit->frac_aligned_query), shift @$d);

521

673

if( $count == 0 ) {

522

674

while( my $hsp = $hit->next_hsp ) {

523

675

ok($hsp->query->start, 1);

527

679

ok($hsp->hit->end, 212);

528

680

ok($hsp->hit->strand, 1);

529

681

ok($hsp->length('hsp'), 60);

530

ok($hsp->evalue == '3e-18');

682

ok(sprintf("%g",$hsp->evalue), sprintf("%g",'3e-18'));

531

683

ok($hsp->score, 48);

532

684

ok($hsp->bits,95.6);

533

685

ok(sprintf("%.2f",$hsp->percent_identity), 96.67);

541

693

ok($hsp->query_string, 'aggaatgctgtttaattggaatcgtacaatggagaatttgacggaaatagaatcaacgat');

542

694

ok($hsp->hit_string, 'aggaatgctgtttaattggaatca-acaatggagaatttgacggaaatagaatcaacgat');

543

695

ok($hsp->homology_string, '||||||||||||||||||||||| |||||||||||||||||||||||||||||||||||');

544

}

545

}

546

last if( $count++ > @valid );

547

}

548

549

# TODO: Flesh this test out!

550

$searchio = new Bio::SearchIO ('-format' => 'psiblast',

696

ok(sprintf("%.2f",$hsp->get_aln->overall_percentage_identity), 96.67);

697

ok(sprintf("%.2f",$hsp->get_aln->percentage_identity), 98.31);

698

}

699

}

700

last if( $count++ > @valid );

701

}

702

703

#WU-BlastX test

704

705

$searchio = new Bio::SearchIO('-format' => 'blast',

706

'-file' => Bio::Root::IO->catfile('t','data','dnaEbsub_ecoli.wublastx'));

707

708

$result = $searchio->next_result;

709

ok($result->database_name, 'ecoli.aa');

710

ok($result->database_letters, 1358990);

711

ok($result->database_entries, 4289);

712

ok($result->algorithm, 'BLASTX');

713

ok($result->algorithm_version, qr/^2\.0MP\-WashU/);

714

ok($result->query_name, 'gi|142864|gb|M10040.1|BACDNAE');

715

ok($result->query_description, 'B.subtilis dnaE gene encoding DNA primase, complete cds');

716

ok($result->query_accession, 'M10040.1');

717

ok($result->query_length, 2001);

718

ok($result->get_parameter('matrix'), 'blosum62');

719

720

ok($result->get_statistic('lambda'), 0.318);

721

ok($result->get_statistic('kappa'), 0.135);

722

ok($result->get_statistic('entropy'),0.401 );

723

724

ok($result->get_statistic('dbentries'), 4289);

725

726

@valid = ( [ 'gi|1789447|gb|AAC76102.1|', 581, 'AAC76102', '1.1e-74', 671]);

727

$count = 0;

728

729

while( my $hit = $result->next_hit ) {

730

my $d = shift @valid;

731

ok($hit->name, shift @$d);

732

ok($hit->length, shift @$d);

733

ok($hit->accession, shift @$d);

734

ok(sprintf("%g",$hit->significance), sprintf("%g",shift @$d) );

735

ok($hit->raw_score, shift @$d );

736

737

if( $count == 0 ) {

738

while( my $hsp = $hit->next_hsp ) {

739

ok($hsp->query->start, 21);

740

ok($hsp->query->end, 1265);

741

ok($hsp->query->strand, 1);

742

ok($hsp->hit->start, 1);

743

ok($hsp->hit->end, 413);

744

ok($hsp->hit->strand, 0);

745

ok($hsp->length('hsp'), 421);

746

ok(sprintf("%g",$hsp->evalue), sprintf("%g",'1.1e-74'));

747

ok(sprintf("%g",$hsp->pvalue), sprintf("%g",'1.1e-74'));

748

ok($hsp->score,671);

749

ok($hsp->bits,265.8);

750

ok(sprintf("%.2f",$hsp->percent_identity), 35.87);

751

ok(sprintf("%.4f",$hsp->frac_identical('query')), 0.1213);

752

ok(sprintf("%.4f",$hsp->frac_identical('hit')), 0.3656);

753

ok(sprintf("%.4f",$hsp->frac_conserved('hsp')), 0.5297);

754

ok($hsp->query->frame(), 2);

755

ok($hsp->hit->frame(), 0);

756

ok($hsp->gaps('query'), 6);

757

ok($hsp->gaps('hit'), 8);

758

ok($hsp->gaps, 14);

759

ok($hsp->query_string, 'MGNRIPDEIVDQVQKSADIVEVIGDYVQLKKQGRNYFGLCPFHGESTPSFSVSPDKQIFHCFGCGAGGNVFSFLRQMEGYSFAESVSHLADKYQIDFPDDITVHSGARP---ESSGEQKMAEAHELLKKFYHHLLINTKEGQEALDYLLSRGFTKELINEFQIGYALDSWDFITKFLVKRGFSEAQMEKAGLLIRREDGSGYFDRFRNRVMFPIHDHHGAVVAFSGRALGSQQPKYMNSPETPLFHKSKLLYNFYKARLHIRKQERAVLFEGFADVYTAVSSDVKESIATMGTSLTDDHVKILRRNVEEIILCYDSDKAGYEATLKASELL---QKKGCKVRVAMIPDGLDPDDYIKKFGGEKFKNDIIDASVTVMAFKMQYFRKGKNLSDEGDRLAYIKDVLKEISTLSGSLEQEVYVKQ');

760

ok($hsp->hit_string, 'MAGRIPRVFINDLLARTDIVDLIDARVKLKKQGKNFHACCPFHNEKTPSFTVNGEKQFYHCFGCGAHGNAIDFLMNYDKLEFVETVEELAAMHNLEVPFE----AGSGPSQIERHQRQTLYQLMDGLNTFYQQSL-QQPVATSARQYLEKRGLSHEVIARFAIGFAPPGWDNVLKRFGGNPENRQSLIDAGMLVTNDQGRSY-DRFRERVMFPIRDKRGRVIGFGGRVLGNDTPKYLNSPETDIFHKGRQLYGLYEAQQDNAEPNRLLVVEGYMDVVALAQYGINYAVASLGTSTTADHIQLLFRATNNVICCYDGDRAGRDAAWRALETALPYMTDGRQLRFMFLPDGEDPDTLVRKEGKEAFEARM-EQAMPLSAFLFNSLMPQVDLSTPDGRARLSTLALPLISQVPGETLR-IYLRQ');

761

ok($hsp->homology_string, 'M RIP ++ + DIV++I V+LKKQG+N+ CPFH E TPSF+V+ +KQ +HCFGCGA GN FL + F E+V LA + ++ P + +G+ P E Q + + + L FY L A YL RG + E+I F IG+A WD + K + + AG+L+ + G Y DRFR RVMFPI D G V+ F GR LG+ PKY+NSPET +FHK + LY Y+A+ + R ++ EG+ DV + ++A++GTS T DH+++L R +I CYD D+AG +A +A E G ++R +PDG DPD ++K G E F+ + + ++ + AF +LS R L IS + G + +Y++Q');

762

}

763

}

764

last if( $count++ > @valid );

765

}

766

767

#WU-TBlastN test

768

769

$searchio = new Bio::SearchIO('-format' => 'blast',

770

'-file' => Bio::Root::IO->catfile('t','data','dnaEbsub_ecoli.wutblastn'));

771

772

$result = $searchio->next_result;

773

ok($result->database_name, 'ecoli.nt');

774

ok($result->database_letters, 4662239);

775

ok($result->database_entries, 400);

776

ok($result->algorithm, 'TBLASTN');

777

ok($result->algorithm_version, qr/^2\.0MP\-WashU/);

778

ok($result->query_name, 'gi|142865|gb|AAA22406.1|');

779

ok($result->query_description, 'DNA primase');

780

ok($result->query_accession, 'AAA22406.1');

781

ok($result->query_length, 603);

782

ok($result->get_parameter('matrix'), 'blosum62');

783

784

ok($result->get_statistic('lambda'), '0.320');

785

ok($result->get_statistic('kappa'), 0.136);

786

ok($result->get_statistic('entropy'),0.387 );

787

788

ok($result->get_statistic('dbentries'), 400);

789

790

@valid = ( [ 'gi|1789441|gb|AE000388.1|AE000388', 10334, 'AE000388', '1.4e-73', 671]);

791

$count = 0;

792

793

while( my $hit = $result->next_hit ) {

794

my $d = shift @valid;

795

ok($hit->name, shift @$d);

796

ok($hit->length, shift @$d);

797

ok($hit->accession, shift @$d);

798

ok(sprintf("%g",$hit->significance), sprintf("%g",shift @$d) );

799

ok($hit->raw_score, shift @$d );

800

801

if( $count == 0 ) {

802

while( my $hsp = $hit->next_hsp ) {

803

ok($hsp->query->start, 1);

804

ok($hsp->query->end, 415);

805

ok($hsp->query->strand, 0);

806

ok($hsp->hit->start, 4778);

807

ok($hsp->hit->end, 6016);

808

ok($hsp->hit->strand, 1);

809

ok($hsp->length('hsp'), 421);

810

ok(sprintf("%g",$hsp->evalue), sprintf("%g",'1.4e-73'));

811

ok(sprintf("%g",$hsp->pvalue), sprintf("%g",'1.4e-73'));

812

ok($hsp->score,671);

813

ok($hsp->bits,265.8);

814

ok(sprintf("%.2f",$hsp->percent_identity), 35.87);

815

ok(sprintf("%.4f",$hsp->frac_identical('hit')), 0.1219);

816

ok(sprintf("%.4f",$hsp->frac_identical('query')), 0.3639);

817

ok(sprintf("%.4f",$hsp->frac_conserved('hsp')), 0.5297);

818

ok($hsp->query->frame(), 0);

819

ok($hsp->hit->frame(), 1);

820

ok($hsp->gaps('query'), 6);

821

ok($hsp->gaps('hit'), 8);

822

ok($hsp->gaps, 14);

823

824

825

826

}

827

}

828

last if( $count++ > @valid );

829

}

830

831

# WU-BLAST TBLASTX

832

$searchio = new Bio::SearchIO('-format' => 'blast',

833

'-file' => Bio::Root::IO->catfile('t','data','dnaEbsub_ecoli.wutblastx'));

834

835

$result = $searchio->next_result;

836

ok($result->database_name, 'ecoli.nt');

837

ok($result->database_letters, 4662239);

838

ok($result->database_entries, 400);

839

ok($result->algorithm, 'TBLASTX');

840

ok($result->algorithm_version, qr/^2\.0MP\-WashU/);

841

ok($result->query_name, 'gi|142864|gb|M10040.1|BACDNAE');

842

ok($result->query_description, 'B.subtilis dnaE gene encoding DNA primase, complete cds');

843

ok($result->query_accession, 'M10040.1');

844

ok($result->query_length, 2001);

845

ok($result->get_parameter('matrix'), 'blosum62');

846

847

ok($result->get_statistic('lambda'), 0.318);

848

ok($result->get_statistic('kappa'), 0.135);

849

ok($result->get_statistic('entropy'),0.401 );

850

ok($result->get_statistic('dbentries'), 400);

851

852

@valid = ( [ 'gi|1789441|gb|AE000388.1|AE000388', 10334, 'AE000388', '6.4e-70', 318],

853

[ 'gi|2367383|gb|AE000509.1|AE000509', 10589, 'AE000509', '0.9992', 59]

854

);

855

$count = 0;

856

857

while( my $hit = $result->next_hit ) {

858

my $d = shift @valid;

859

ok($hit->name, shift @$d);

860

ok($hit->length, shift @$d);

861

ok($hit->accession, shift @$d);

862

# using e here to deal with 0.9992 coming out right here as well

863

ok(sprintf("%g",$hit->significance), sprintf("%g",shift @$d) );

864

ok($hit->raw_score, shift @$d );

865

866

if( $count == 0 ) {

867

my $hspcounter = 0;

868

while( my $hsp = $hit->next_hsp ) {

869

$hspcounter++;

870

if( $hspcounter == 3 ) {

871

# let's actually look at the 3rd HSP

872

ok($hsp->query->start, 441);

873

ok($hsp->query->end, 617);

874

ok($hsp->query->strand, 1);

875

ok($hsp->hit->start, 5192);

876

ok($hsp->hit->end, 5368);

877

ok($hsp->hit->strand, 1);

878

ok($hsp->length('hsp'), 59);

879

ok(sprintf("%g",$hsp->evalue), sprintf("%g",'6.4e-70'));

880

ok(sprintf("%g",$hsp->pvalue), sprintf("%g",'6.4e-70'));

881

ok($hsp->score,85);

882

ok($hsp->bits,41.8);

883

ok(sprintf("%.2f",$hsp->percent_identity), '32.20');

884

ok(sprintf("%.4f",$hsp->frac_identical('hit')), 0.1073);

885

ok(sprintf("%.4f",$hsp->frac_identical('query')), 0.1073);

886

ok(sprintf("%.4f",$hsp->frac_conserved('hsp')), 0.4746);

887

ok($hsp->query->frame(), 2);

888

ok($hsp->hit->frame(), 1);

889

ok($hsp->gaps('query'), 0);

890

ok($hsp->gaps('hit'), 0);

891

ok($hsp->gaps, 0);

892

ok($hsp->query_string, 'ALDYLLSRGFTKELINEFQIGYALDSWDFITKFLVKRGFSEAQMEKAGLLIRREDGSGY');

893

ok($hsp->hit_string, 'ARQYLEKRGLSHEVIARFAIGFAPPGWDNVLKRFGGNPENRQSLIDAGMLVTNDQGRSY');

894

ok($hsp->homology_string, 'A YL RG + E+I F IG+A WD + K + + AG+L+ + G Y');

895

}

896

}

897

} elsif( $count == 1 ) {

898

while( my $hsp = $hit->next_hsp ) {

899

ok($hsp->query->start, 587);

900

ok($hsp->query->end, 706);

901

ok($hsp->query->strand, -1);

902

ok($hsp->hit->start, 4108);

903

ok($hsp->hit->end, 4227);

904

ok($hsp->hit->strand, -1);

905

ok($hsp->length('hsp'), 40);

906

ok($hsp->evalue == '7.1');

907

ok($hsp->pvalue == '1.00');

908

ok($hsp->score,59);

909

ok($hsp->bits,29.9);

910

ok(sprintf("%.2f",$hsp->percent_identity), '37.50');

911

ok(sprintf("%.4f",$hsp->frac_identical('hit')), '0.1250');

912

ok(sprintf("%.4f",$hsp->frac_identical('query')), '0.1250');

913

ok(sprintf("%.4f",$hsp->frac_conserved('hsp')), '0.4750');

914

ok($hsp->query->frame(), 2);

915

ok($hsp->hit->frame(), 2);

916

ok($hsp->gaps('query'), 0);

917

ok($hsp->gaps('hit'), 0);

918

ok($hsp->gaps, 0);

919

ok($hsp->query_string, 'WLPRALPEKATTAP**SWIGNMTRFLKRSKYPLPSSRLIR');

920

ok($hsp->hit_string, 'WLSRTTVGSSTVSPRTFWITRMKVKLSSSKVTLPSTKSTR');

921

ok($hsp->homology_string, 'WL R +T +P WI M L SK LPS++ R');

922

last;

923

}

924

}

925

last if( $count++ > @valid );

926

}

927

928

# Do a multiblast report test

929

$searchio = new Bio::SearchIO ('-format' => 'blast',

930

'-file' => Bio::Root::IO->catfile('t','data','multi_blast.bls'));

931

932

my @expected = qw(CATH_RAT CATL_HUMAN CATL_RAT PAPA_CARPA);

933

while( my $result = $searchio->next_result ) {

934

ok($result->query_name, shift @expected, "Multiblast query test");

935

}

936

937

938

# Test GCGBlast parsing

939

940

$searchio = new Bio::SearchIO('-format' => 'blast',

941

'-file' => Bio::Root::IO->catfile('t','data', 'test.gcgblast'));

942

$result = $searchio->next_result();

943

944

ok($result->query_name, '/v0/people/staji002/test.gcg');

945

ok($result->algorithm, 'BLASTP');

946

ok($result->algorithm_version, '2.2.1 [Apr-13-2001]');

947

ok($result->database_name, 'pir');

948

ok($result->database_entries, 274514);

949

ok($result->database_letters, 93460074);

950

951

$hit = $result->next_hit;

952

ok($hit->name, 'PIR2:S44629');

953

ok($hit->length, 628);

954

ok($hit->accession, 'PIR2:S44629');

955

skip('Significance parsing broken for GCG-BLAST Hits -- see HSP',$hit->significance, '2e-08' );

956

skip('Raw score parsing broken for GCG-BLAST Hits -- see HSP',$hit->raw_score, 57 );

957

958

$hsp = $hit->next_hsp;

959

ok(sprintf("%g",$hsp->evalue), sprintf("%g",'2e-08'));

960

ok($hsp->bits, '57.0');

961

ok($hsp->score, 136);

962

ok(int($hsp->percent_identity), 28);

963

ok(sprintf("%.2f",$hsp->frac_identical('query')), 0.29);

964

ok($hsp->frac_conserved('total'), 69/135);

965

ok($hsp->gaps('total'), 8);

966

ok($hsp->gaps('hit'), 6);

967

ok($hsp->gaps('query'), 2);

968

969

ok($hsp->hit->start, 342);

970

ok($hsp->hit->end, 470);

971

ok($hsp->query->start, 3);

972

ok($hsp->query->end, 135);

973

974

ok($hsp->query_string, 'CAAEFDFMEKETPLRYTKTXXXXXXXXXXXXXXRKIISDMWGVLAKQQTHVRKHQFDHGELVYHALQLLAYTALGILIMRLKLFLTPYMCVMASLICSRQLFGW--LFCKVHPGAIVFVILAAMSIQGSANLQTQ');

975

ok($hsp->hit_string, 'CSAEFDFIQYSTIEKLCGTLLIPLALISLVTFVFNFVKNT-NLLWRNSEEIG----ENGEILYNVVQLCCSTVMAFLIMRLKLFMTPHLCIVAALFANSKLLGGDRISKTIRVSALVGVI-AILFYRGIPNIRQQ');

976

ok($hsp->homology_string, 'C+AEFDF++ T + T + + +L + + ++GE++Y+ +QL T + LIMRLKLF+TP++C++A+L + +L G + + A+V VI A + +G N++ Q');

977

978

979

$searchio = new Bio::SearchIO ('-format' => 'blast',

551

980

'-file' => Bio::Root::IO->catfile('t','data','HUMBETGLOA.tblastx'));

552

981

553

982

$result = $searchio->next_result;

554

983

555

984

ok($result);

556

985

$hit = $result->next_hit;

986

ok($hit->accession, 'AE000479');

987

ok($hit->bits, 33.6);

557

988

$hsp = $hit->next_hsp;

989

ok($hit->hsp->bits,$hsp->bits);

990

558

991

ok($hsp->get_aln->isa('Bio::Align::AlignI'));

559

992

my $writer = Bio::SearchIO::Writer::HitTableWriter->new(

560

993

-columns => [qw(

562

995

query_length

563

996

hit_name

564

997

hit_length

998

bits

999

score

565

1000

frac_identical_query

566

1001

expect

567

1002

)] );

568

1003

569

1004

my $out = new Bio::SearchIO(-writer => $writer,

570

-file => ">searchio.out");

1005

-file => ">searchio.out");

571

1006

$out->write_result($result, 1);

572

1007

ok(-e 'searchio.out');

573

1008

my $writerhtml = new Bio::SearchIO::Writer::HTMLResultWriter();

574

1009

my $outhtml = new Bio::SearchIO(-writer => $writerhtml,

575

1010

-file => ">searchio.html");

1011

$outhtml->write_result($result, 1);

576

1012

ok(-e "searchio.html");

577

1013

578

END {

579

unlink 'searchio.out';

580

unlink 'searchio.html';

581

582

}

1014

unlink 'searchio.out';

1015

unlink 'searchio.html';

1016

1017

#test all the database accession number formats

1018

$searchio = new Bio::SearchIO(-format => 'blast',

1019

-file => 't/data/testdbaccnums.out');

1020

$result = $searchio->next_result;

1021

1022

@valid = ( ['pir||T14789','T14789','T14789','CAB53709','AAH01726'],

1023

['gb|NP_065733.1|CYT19', 'NP_065733','CYT19'],

1024

['emb|XP_053690.4|Cyt19','XP_053690'],

1025

['dbj|NP_056277.2|DKFZP586L0724','NP_056277'],

1026

['prf||XP_064862.2','XP_064862'],

1027

['pdb|BAB13968.1|1','BAB13968'],

1028

['sp|Q16478|GLK5_HUMAN','Q16478'],

1029

['pat|US|NP_002079.2','NP_002079'],

1030

['bbs|NP_079463.2|','NP_079463'],

1031

['gnl|db1|NP_002444.1','NP_002444'],

1032

['ref|XP_051877.1|','XP_051877'],

1033

['lcl|AAH16829.1|','AAH16829'],

1034

['gi|1|gb|NP_065733.1|CYT19','NP_065733'],

1035

['gi|2|emb|XP_053690.4|Cyt19','XP_053690'],

1036

['gi|3|dbj|NP_056277.2|DKFZP586L0724','NP_056277'],

1037

['gi|4|pir||T14789','T14789'],

1038

['gi|5|prf||XP_064862.2','XP_064862'],

1039

['gi|6|pdb|BAB13968.1|1','BAB13968'],

1040

['gi|7|sp|Q16478|GLK5_HUMAN','Q16478'],

1041

['gi|8|pat|US|NP_002079.2','NP_002079'],

1042

['gi|9|bbs|NP_079463.2|','NP_079463'],

1043

['gi|10|gnl|db1|NP_002444.1','NP_002444'],

1044

['gi|11|ref|XP_051877.1|','XP_051877'],

1045

['gi|12|lcl|AAH16829.1|','AAH16829'],

1046

['MY_test_ID','MY_test_ID']

1047

);

1048

1049

$hit = $result->next_hit;

1050

my $d = shift @valid;

1051

ok($hit->name, shift @$d);

1052

ok($hit->accession, shift @$d);

1053

my @accnums = $hit->each_accession_number;

1054

foreach my $a (@accnums) {

1055

ok($a, shift @$d);

1056

}

1057

$d = shift @valid;

1058

$hit = $result->next_hit;

1059

ok($hit->name, shift @$d);

1060

ok($hit->accession, shift @$d);

1061

ok($hit->locus, shift @$d);

1062

1063

while( $hit = $result->next_hit ) {

1064

my $d = shift @valid;

1065

ok($hit->name, shift @$d);

1066

ok($hit->accession, shift @$d);

1067

}

1068

1069

# Parse MEGABLAST

1070

1071

# parse the BLAST-like output

1072

my $infile = Bio::Root::IO->catfile(qw(t data 503384.MEGABLAST.2));

1073

my $in = new Bio::SearchIO(-file => $infile,

1074

-format => 'blast'); # this is megablast

1075

# blast-like output

1076

my $r = $in->next_result;

1077

my @dcompare = ( ['Contig3700', 5631, 785, '0.0', 785, '0.0', 396, 639, 12,

1078

8723,9434, 1, 4083, 4794, -1],

1079

['Contig3997', 12734, 664, '0.0', 664, '0.0', 335, 401, 0,

1080

1282, 1704, 1, 1546, 1968,-1 ],

1081

['Contig634', 858, 486, '1e-136', 486, '1e-136', 245, 304, 3,

1082

7620, 7941, 1, 1, 321, -1],

1083

['Contig1853', 2314, 339, '1e-91',339, '1e-91', 171, 204, 0,

1084

6406, 6620, 1, 1691, 1905, 1]

1085

);

1086

1087

ok($r->query_name, '503384');

1088

ok($r->query_description, '11337 bp 2 contigs');

1089

ok($r->query_length, 11337);

1090

ok($r->database_name, 'cneoA.nt ');

1091

ok($r->database_letters, 17206226);

1092

ok($r->database_entries, 4935);

1093

1094

while( my $hit = $r->next_hit ) {

1095

my $d = shift @dcompare;

1096

ok($hit->name, shift @$d);

1097

ok($hit->length, shift @$d);

1098

ok($hit->raw_score, shift @$d);

1099

ok($hit->significance, shift @$d);

1100

1101

my $hsp = $hit->next_hsp;

1102

ok($hsp->bits, shift @$d);

1103

ok($hsp->evalue, shift @$d);

1104

ok($hsp->score, shift @$d);

1105

ok($hsp->num_identical, shift @$d);

1106

ok($hsp->gaps('total'), shift @$d);

1107

ok($hsp->query->start, shift @$d);

1108

ok($hsp->query->end, shift @$d);

1109

ok($hsp->query->strand, shift @$d);

1110

ok($hsp->hit->start, shift @$d);

1111

ok($hsp->hit->end, shift @$d);

1112

ok($hsp->hit->strand, shift @$d);

1113

}

1114

1115

1116

# parse the another megablast format

1117

1118

$infile = Bio::Root::IO->catfile(qw(t data 503384.MEGABLAST.0));

1119

1120

# this is megablast output type 0

1121

$in = new Bio::SearchIO(-file => $infile,

1122

-report_format => 0,

1123

-format => 'megablast');

1124

$r = $in->next_result;

1125

@dcompare = (

1126

['Contig634', 7620, 7941, 1, 1, 321, -1],

1127

['Contig1853', 6406, 6620, 1, 1691, 1905, 1],

1128

['Contig3700', 8723,9434, 1, 4083, 4794, -1],

1129

['Contig3997', 1282, 1704, 1, 1546, 1968,-1 ],

1130

);

1131

1132

ok($r->query_name, '503384');

1133

1134

while( my $hit = $r->next_hit ) {

1135

my $d = shift @dcompare;

1136

ok($hit->name, shift @$d);

1137

my $hsp = $hit->next_hsp;

1138

ok($hsp->query->start, shift @$d);

1139

ok($hsp->query->end, shift @$d);

1140

ok($hsp->query->strand, shift @$d);

1141

ok($hsp->hit->start, shift @$d);

1142

ok($hsp->hit->end, shift @$d);

1143

ok($hsp->hit->strand, shift @$d);

1144

}

1145

1146

1147

1148

# Let's test RPS-BLAST

1149

1150

my $parser = new Bio::SearchIO(-format => 'blast',

1151

-file => Bio::Root::IO->catfile(qw(t data ecoli_domains.rpsblast)));

1152

1153

$r = $parser->next_result;

1154

ok($r->query_name, 'gi|1786183|gb|AAC73113.1|');

1155

ok($r->num_hits, 7);

1156

$hit = $r->next_hit;

1157

ok($hit->name, 'gnl|CDD|3919');

1158

ok($hit->significance, 0.064);

1159

ok($hit->raw_score, 28);

1160

$hsp = $hit->next_hsp;

1161

ok($hsp->query->start, 599);

1162

ok($hsp->query->end,655);

1163

ok($hsp->hit->start,23);

1164

ok($hsp->hit->end,76);

1165

1166

1167

# Test PSI-BLAST parsing

1168

1169

$searchio = new Bio::SearchIO ('-format' => 'blast',

1170

'-file' => Bio::Root::IO->catfile('t','data','psiblastreport.out'));

1171

1172

$result = $searchio->next_result;

1173

1174

ok($result->database_name, '/home/peter/blast/data/swissprot.pr');

1175

ok($result->database_entries, 88780);

1176

ok($result->database_letters, 31984247);

1177

1178

ok($result->algorithm, 'BLASTP');

1179

ok($result->algorithm_version, qr/^2\.0\.14/);

1180

ok($result->query_name, 'CYS1_DICDI');

1181

ok($result->query_length, 343);

1182

ok($result->get_statistic('kappa') == 0.0491);

1183

ok($result->get_statistic('lambda') == 0.270);

1184

ok($result->get_statistic('entropy') == 0.230);

1185

ok($result->get_statistic('dbletters'), 31984247);

1186

ok($result->get_statistic('dbentries'), 88780);

1187

ok($result->get_statistic('effective_hsplength'), 49);

1188

ok($result->get_statistic('effectivespace'), 8124403938);

1189

ok($result->get_parameter('matrix'), 'BLOSUM62');

1190

ok($result->get_parameter('gapopen'), 11);

1191

ok($result->get_parameter('gapext'), 1);

1192

1193

my @valid_hit_data = ( [ 'sp|P04988|CYS1_DICDI', 343, 'P04988', '0', 721],

1194

[ 'sp|P43295|A494_ARATH', 313, 'P43295', '1e-75', 281],

1195

[ 'sp|P25804|CYSP_PEA', 363, 'P25804', '1e-74', 278]);

1196

my @valid_iter_data = ( [ 127, 127, 0, 109, 18, 0, 0, 0, 0],

1197

[ 157, 40, 117, 2, 38, 0, 109, 3, 5]);

1198

my $iter_count = 0;

1199

while( $iter = $result->next_iteration ) {

1200

$iter_count++;

1201

my $di = shift @valid_iter_data;

1202

ok($iter->number, $iter_count);

1203

1204

ok($iter->num_hits, shift @$di);

1205

ok($iter->num_hits_new, shift @$di);

1206

ok($iter->num_hits_old, shift @$di);

1207

ok(scalar($iter->newhits_below_threshold), shift @$di);

1208

ok(scalar($iter->newhits_not_below_threshold), shift @$di);

1209

ok(scalar($iter->newhits_unclassified), shift @$di);

1210

ok(scalar($iter->oldhits_below_threshold), shift @$di);

1211

ok(scalar($iter->oldhits_newly_below_threshold), shift @$di);

1212

ok(scalar($iter->oldhits_not_below_threshold), shift @$di);

1213

1214

my $hit_count = 0;

1215

if ($iter_count == 1) {

1216

while( $hit = $result->next_hit ) {

1217

my $d = shift @valid_hit_data;

1218

1219

ok($hit->name, shift @$d);

1220

ok($hit->length, shift @$d);

1221

ok($hit->accession, shift @$d);

1222

ok(sprintf("%g",$hit->significance), sprintf("%g",shift @$d) );

1223

ok($hit->bits, shift @$d );

1224

1225

if( $hit_count == 1 ) {

1226

while( my $hsp = $hit->next_hsp ){

1227

ok($hsp->query->start, 32);

1228

ok($hsp->query->end, 340);

1229

ok($hsp->hit->start, 3);

1230

ok($hsp->hit->end, 307);

1231

ok($hsp->length('hsp'), 316);

1232

ok($hsp->start('hit'), $hsp->hit->start);

1233

ok($hsp->end('query'), $hsp->query->end);

1234

ok($hsp->strand('sbjct'), $hsp->subject->strand);# alias for hit

1235

ok($hsp->evalue == 1e-75);

1236

ok($hsp->score, 712);

1237

ok($hsp->bits, 281);

1238

ok(sprintf("%.1f",$hsp->percent_identity), 46.5);

1239

ok(sprintf("%.4f",$hsp->frac_identical('query')), 0.4757);

1240

ok(sprintf("%.3f",$hsp->frac_identical('hit')), 0.482);

1241

ok($hsp->gaps, 18);

1242

}

1243

}

1244

last if( $hit_count++ > @valid_hit_data );

1245

}

1246

}

1247

}

1248

1249

# Test filtering

1250

1251

$searchio = new Bio::SearchIO ( '-format' => 'blast',

1252

'-file' => Bio::Root::IO->catfile('t','data','ecolitst.bls'),

1253

'-signif' => 1e-100);

1254

1255

@valid = qw(gb|AAC73113.1|);

1256

$r = $searchio->next_result;

1257

1258

while( my $hit = $r->next_hit ) {

1259

ok($hit->name, shift @valid);

1260

}

1261

1262

$searchio = new Bio::SearchIO ( '-format' => 'blast',

1263

'-file' => Bio::Root::IO->catfile('t','data','ecolitst.bls'),

1264

'-score' => 100);

1265

1266

@valid = qw(gb|AAC73113.1| gb|AAC76922.1| gb|AAC76994.1|);

1267

$r = $searchio->next_result;

1268

1269

while( my $hit = $r->next_hit ) {

1270

ok($hit->name, shift @valid);

1271

}

1272

1273

$searchio = new Bio::SearchIO ( '-format' => 'blast',

1274

'-file' => Bio::Root::IO->catfile('t','data','ecolitst.bls'),

1275

'-bits' => 200);

1276

1277

@valid = qw(gb|AAC73113.1| gb|AAC76922.1|);

1278

$r = $searchio->next_result;

1279

1280

while( my $hit = $r->next_hit ) {

1281

ok($hit->name, shift @valid);

1282

}

1283

1284

1285

my $filt_func = sub{ my $hit=shift;

1286

$hit->frac_identical('query') >= 0.31

1287

};

1288

1289

$searchio = new Bio::SearchIO ( '-format' => 'blast',

1290

'-file' => Bio::Root::IO->catfile('t','data','ecolitst.bls'),

1291

'-hit_filter' => $filt_func);

1292

1293

@valid = qw(gb|AAC73113.1| gb|AAC76994.1|);

1294

$r = $searchio->next_result;

1295

1296

while( my $hit = $r->next_hit ) {

1297

ok($hit->name, shift @valid);

1298

}

1299

1300

1301

1302

1303

# bl2seq parsing testing

1304

1305

# this is blastp bl2seq

1306

$searchio = new Bio::SearchIO(-format => 'blast',

1307

-file => Bio::Root::IO->catfile(qw(t data

1308

bl2seq.out)));

1309

$result = $searchio->next_result;

1310

ok($result);

1311

ok($result->query_name, '');

1312

ok($result->algorithm, 'BLASTP');

1313

$hit = $result->next_hit;

1314

ok($hit->name, 'ALEU_HORVU');

1315

ok($hit->length, 362);

1316

$hsp = $hit->next_hsp;

1317

ok($hsp->score, 481);

1318

ok($hsp->bits, 191);

1319

ok(int $hsp->percent_identity, 34);

1320

ok($hsp->evalue, '2e-53');

1321

ok(int($hsp->frac_conserved*$hsp->length), 167);

1322

ok($hsp->query->start, 28);

1323

ok($hsp->query->end, 343);

1324

ok($hsp->hit->start, 60);

1325

ok($hsp->hit->end,360);

1326

ok($hsp->gaps, 27);

1327

1328

# this is blastn bl2seq

1329

$searchio = new Bio::SearchIO(-format => 'blast',

1330

-file => Bio::Root::IO->catfile

1331

(qw(t data bl2seq.blastn.rev)));

1332

$result = $searchio->next_result;

1333

ok($result);

1334

ok($result->query_name, '');

1335

ok($result->algorithm, 'BLASTN');

1336

ok($result->query_length, 180);

1337

$hit = $result->next_hit;

1338

ok($hit->length, 179);

1339

ok($hit->name, 'human');

1340

$hsp = $hit->next_hsp;

1341

ok($hsp->score, 27);

1342

ok($hsp->bits, '54.0');

1343

ok(int $hsp->percent_identity, 88);

1344

ok($hsp->evalue, '2e-12');

1345

ok(int($hsp->frac_conserved*$hsp->length), 83);

1346

ok($hsp->query->start, 94);

1347

ok($hsp->query->end, 180);

1348

ok($hsp->query->strand, 1);

1349

ok($hsp->hit->strand, -1);

1350

ok($hsp->hit->start, 1);

1351

ok($hsp->hit->end,94);

1352

ok($hsp->gaps, 7);

1353

1354

# this is blastn bl2seq

1355

$searchio = new Bio::SearchIO(-format => 'blast',

1356

-file => Bio::Root::IO->catfile

1357

(qw(t data bl2seq.blastn)));

1358

$result = $searchio->next_result;

1359

ok($result);

1360

ok($result->query_name, '');

1361

ok($result->query_length, 180);

1362

ok($result->algorithm, 'BLASTN');

1363

$hit = $result->next_hit;

1364

ok($hit->name, 'human');

1365

ok($hit->length, 179);

1366

$hsp = $hit->next_hsp;

1367

ok($hsp->score, 27);

1368

ok($hsp->bits, '54.0');

1369

ok(int $hsp->percent_identity, 88);

1370

ok($hsp->evalue, '2e-12');

1371

ok(int($hsp->frac_conserved*$hsp->length), 83);

1372

ok($hsp->query->start, 94);

1373

ok($hsp->query->end, 180);

1374

ok($hsp->query->strand, 1);

1375

ok($hsp->hit->strand, 1);

1376

ok($hsp->hit->start, 86);

1377

ok($hsp->hit->end,179);

1378

ok($hsp->gaps, 7);

1379

1380

# this is blastp bl2seq

1381

$searchio = new Bio::SearchIO(-format => 'blast',

1382

-file => Bio::Root::IO->catfile

1383

(qw(t data bl2seq.bug940.out)));

1384

$result = $searchio->next_result;

1385

ok($result);

1386

ok($result->query_name, 'zinc');

1387

ok($result->algorithm, 'BLASTP');

1388

ok($result->query_description, 'finger protein 135 (clone pHZ-17) [Homo sapiens]. neo_id RS.ctg14243-000000.6.0');

1389

ok($result->query_length, 469);

1390

$hit = $result->next_hit;

1391

ok($hit->name, 'gi|4507985|');

1392

ok($hit->description,'zinc finger protein 135 (clone pHZ-17) [Homo sapiens]. neo_id RS.ctg14243-000000.6.0');

1393

ok($hit->length, 469);

1394

$hsp = $hit->next_hsp;

1395

ok($hsp->score, 1626);

1396

ok($hsp->bits, 637);

1397

ok(int $hsp->percent_identity, 66);

1398

ok($hsp->evalue, '0.0');

1399

ok(int($hsp->frac_conserved*$hsp->length), 330);

1400

ok($hsp->query->start, 121);

1401

ok($hsp->query->end, 469);

1402

ok($hsp->hit->start, 1);

1403

ok($hsp->hit->end,469);

1404

ok($hsp->gaps, 120);

1405

ok($hit->next_hsp); # there is more than one HSP here,

1406

# make sure it is parsed at least

1407

1408

# cannot distinguish between blastx and tblastn reports

1409

# so we're only testing a blastx report for now

1410

1411

# this is blastx bl2seq

1412

$searchio = new Bio::SearchIO(-format => 'blast',

1413

-file => Bio::Root::IO->catfile

1414

(qw(t data bl2seq.blastx.out)));

1415

$result = $searchio->next_result;

1416

ok($result);

1417

ok($result->query_name, 'AE000111.1');

1418

ok($result->query_description, 'Escherichia coli K-12 MG1655 section 1 of 400 of the complete genome');

1419

ok($result->algorithm, 'BLASTX');

1420

ok($result->query_length, 720);

1421

$hit = $result->next_hit;

1422

ok($hit->name, 'AK1H_ECOLI');

1423

ok($hit->description,'P00561 Bifunctional aspartokinase/homoserine dehydrogenase I (AKI-HDI) [Includes: Aspartokinase I ; Homoserine dehydrogenase I ]');

1424

ok($hit->length, 820);

1425

$hsp = $hit->next_hsp;

1426

ok($hsp->score, 634);

1427

ok($hsp->bits, 248);

1428

ok(int $hsp->percent_identity, 100);

1429

ok($hsp->evalue, '2e-70');

1430

ok(int($hsp->frac_conserved*$hsp->length), 128);

1431

ok($hsp->query->start, 1);

1432

ok($hsp->query->end, 384);

1433

ok($hsp->hit->start, 1);

1434

ok($hsp->hit->end,128);

1435

ok($hsp->gaps, 0);

1436

ok($hsp->query->frame,0);

1437

ok($hsp->hit->frame,0);

1438

ok($hsp->query->strand,-1);

1439

ok($hsp->hit->strand,0);

1440

1441

# this is tblastx bl2seq (self against self)

1442

$searchio = new Bio::SearchIO(-format => 'blast',

1443

-file => Bio::Root::IO->catfile

1444

(qw(t data bl2seq.tblastx.out)));

1445

$result = $searchio->next_result;

1446

ok($result);

1447

ok($result->query_name, 'Escherichia');

1448

ok($result->algorithm, 'TBLASTX');

1449

ok($result->query_description, 'coli K-12 MG1655 section 1 of 400 of the complete genome');

1450

ok($result->query_length, 720);

1451

$hit = $result->next_hit;

1452

ok($hit->name, 'gi|1786181|gb|AE000111.1|AE000111');

1453

1454

ok($hit->description,'Escherichia coli K-12 MG1655 section 1 of 400 of the complete genome');

1455

ok($hit->length, 720);

1456

$hsp = $hit->next_hsp;

1457

ok($hsp->score, 1118);

1458

ok($hsp->bits, 515);

1459

ok(int $hsp->percent_identity, 95);

1460

ok($hsp->evalue, '1e-151');

1461

ok(int($hsp->frac_conserved*$hsp->length), 229);

1462

ok($hsp->query->start, 1);

1463

ok($hsp->query->end, 720);

1464

ok($hsp->hit->start, 1);

1465

ok($hsp->hit->end,720);

1466

ok($hsp->gaps, 0);

1467

ok($hsp->query->frame,0);

1468

ok($hsp->hit->frame,0);

1469

ok($hsp->query->strand,1);

1470

ok($hsp->hit->strand,1);

1471

1472

# test blasttable output

1473

my @eqset = qw(

1474

1475

c200-vs-yeast.BLASTN.m9);

1476

$searchio = new Bio::SearchIO(-file => Bio::Root::IO->catfile

1477

(qw(t data c200-vs-yeast.BLASTN)),

1478

-format => 'blast');

1479

$result = $searchio->next_result;

1480

ok($result);

1481

my %ref = &result2hash($result);

1482

ok( scalar keys %ref, 67);

1483

$searchio = new Bio::SearchIO(-file => Bio::Root::IO->catfile

1484

(qw(t data c200-vs-yeast.BLASTN.m8)),

1485

-program_name => 'BLASTN',

1486

-format => 'blasttable');

1487

$result = $searchio->next_result;

1488

my %tester = &result2hash($result);

1489

foreach my $key ( sort keys %ref ) {

1490

ok($tester{$key}, $ref{$key},$key);

1491

}

1492

1493

1494

1495

1496

# Test Blast parsing with B=0 (WU-BLAST)

1497

$searchio = new Bio::SearchIO(-file => Bio::Root::IO->catfile

1498

(qw(t data no_hsps.blastp)),

1499

-format => 'blast');

1500

$result = $searchio->next_result;

1501

ok($result->query_name, 'mgri:MG00189.3');

1502

$hit = $result->next_hit;

1503

ok($hit->name, 'mgri:MG00189.3');

1504

ok($hit->description, 'hypothetical protein 6892 8867 +');

1505

ok($hit->score, 3098);

1506

ok($hit->significance, '0.');

1507

1508

$hit = $result->next_hit;

1509

ok($hit->name, 'fgram:FG01141.1');

1510

ok($hit->description, 'hypothetical protein 47007 48803 -');

1511

ok($hit->score, 2182);

1512

ok($hit->significance, '4.2e-226');

1513

ok($result->num_hits, 415);

1514

# Let's now test if _guess_format is doing its job correctly

1515

my %pair = ( 'filename.blast' => 'blast',

1516

'filename.bls' => 'blast',

1517

'f.blx' => 'blast',

1518

'f.tblx' => 'blast',

1519

'fast.bls' => 'blast',

1520

'f.fasta' => 'fasta',

1521

'f.fa' => 'fasta',

1522

'f.fx' => 'fasta',

1523

'f.fy' => 'fasta',

1524

'f.ssearch' => 'fasta',

1525

'f.SSEARCH.m9' => 'fasta',

1526

'f.m9' => 'fasta',

1527

'f.psearch' => 'fasta',

1528

'f.osearch' => 'fasta',

1529

'f.exon' => 'exonerate',

1530

'f.exonerate' => 'exonerate',

1531

'f.blastxml' => 'blastxml',

1532

'f.xml' => 'blastxml');

1533

while( my ($file,$expformat) = each %pair ) {

1534

ok(Bio::SearchIO->_guess_format($file),$expformat, "$expformat for $file");

1535

}

1536

1537

1538

# Test Wes Barris's reported bug when parsing blastcl3 output which

1539

# has integer overflow

1540

1541

$searchio = new Bio::SearchIO(-file => Bio::Root::IO->catfile

1542

(qw(t data hsinsulin.blastcl3.blastn)),

1543

-format => 'blast');

1544

$result = $searchio->next_result;

1545

ok($result->query_name, 'human');

1546

ok($result->database_letters(), '-24016349');

1547

# this is of course not the right length, but is the what blastcl3

1548

# reports, the correct value is

1549

ok($result->get_statistic('dbletters'),'192913178');

1550

ok($result->get_statistic('dbentries'),'1867771');

1551

1552

1553

# some utilities

1554

# a utility function for comparing result objects

1555

sub result2hash {

1556

my ($result) = @_;

1557

my %hash;

1558

$hash{'query_name'} = $result->query_name;

1559

my $hitcount = 1;

1560

my $hspcount = 1;

1561

foreach my $hit ( $result->hits ) {

1562

$hash{"hit$hitcount\_name"} = $hit->name;

1563

# only going to test order of magnitude

1564

# too hard as these don't always match

1565

# $hash{"hit$hitcount\_signif"} =

1566

# ( sprintf("%.0e",$hit->significance) =~ /e\-?(\d+)/ );

1567

$hash{"hit$hitcount\_bits"} = sprintf("%d",$hit->bits);

1568

1569

foreach my $hsp ( $hit->hsps ) {

1570

$hash{"hsp$hspcount\_bits"} = sprintf("%d",$hsp->bits);

1571

# only going to test order of magnitude

1572

# too hard as these don't always match

1573

# $hash{"hsp$hspcount\_evalue"} =

1574

# ( sprintf("%.0e",$hsp->evalue) =~ /e\-?(\d+)/ );

1575

$hash{"hsp$hspcount\_qs"} = $hsp->query->start;

1576

$hash{"hsp$hspcount\_qe"} = $hsp->query->end;

1577

$hash{"hsp$hspcount\_qstr"} = $hsp->query->strand;

1578

$hash{"hsp$hspcount\_hs"} = $hsp->hit->start;

1579

$hash{"hsp$hspcount\_he"} = $hsp->hit->end;

1580

$hash{"hsp$hspcount\_hstr"} = $hsp->hit->strand;

1581

1582

#$hash{"hsp$hspcount\_pid"} = sprintf("%d",$hsp->percent_identity);

1583

#$hash{"hsp$hspcount\_fid"} = sprintf("%.2f",$hsp->frac_identical);

1584

$hash{"hsp$hspcount\_gaps"} = $hsp->gaps('total');

1585

$hspcount++;

1586

}

1587

$hitcount++;

1588

}

1589

return %hash;

1590

}

1591

1592

1593

__END__

1594

1595

Useful for debugging:

1596

1597

if ($iter_count == 3) {

1598

print "NEWHITS:\n";

1599

foreach ($iter->newhits) {

1600

print " " . $_->name . "\n";

1601

}

1602

print "\nOLDHITS:\n";

1603

foreach ($iter->oldhits) {

1604

print " " . $_->name . "\n";

1605

}

1606

print "\nNEWHITS BELOW:\n";

1607

foreach ($iter->newhits_below_threshold) {

1608

print " " . $_->name . "\n";

1609

}

1610

print "\nNEWHITS NOT BELOW:\n";

1611

foreach ($iter->newhits_not_below_threshold) {

1612

print " " . $_->name . "\n";

1613

}

1614

print "\nNEWHITS UNCLASSIFIED:\n";

1615

foreach ($iter->newhits_unclassified) {

1616

print " " . $_->name . "\n";

1617

}

1618

print "\nOLDHITS BELOW:\n";

1619

foreach ($iter->oldhits_below_threshold) {

1620

print " " . $_->name . "\n";

1621

}

1622

print "\nOLDHITS NEWLY BELOW:\n";

1623

foreach ($iter->oldhits_newly_below_threshold) {

1624

print " " . $_->name . "\n";

1625

}

1626

print "\nOLDHITS NOT BELOW:\n";

1627

foreach ($iter->oldhits_not_below_threshold) {

1628

print " " . $_->name . "\n";

1629

}

1630

}

1631

1632

Older »