Repository: acorg/dark-matter

Issue:
Add read filtering by percent identity match (in DIAMOND output)

Discussion:

Patch:
diff --git a/CHANGELOG.md b/CHANGELOG.md
index d82887e..cbfda73 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -1,3 +1,10 @@
+## 3.1.22 Dec 3, 2019
+
+Added `--percentageIdenticalCutoff` argument to `alignment-panel-civ.py`
+and `noninteractive-alignment-panel.py` and all that implies, down to
+reading DIAMOND output with the `pident` value, storing it, restoring it,
+and filtering on it in `dark.alignments.ReadsAlignmentsFilter`, plus tests.
+
 ## 3.1.21 Nov 30, 2019
 
 Minor changes to HTML output.
diff --git a/bin/alignment-panel-civ.py b/bin/alignment-panel-civ.py
index 03c66c3..decc2c9 100755
--- a/bin/alignment-panel-civ.py
+++ b/bin/alignment-panel-civ.py
@@ -115,12 +115,12 @@ if __name__ == '__main__':
 
     # Args for filtering on ReadsAlignments.
     parser.add_argument(
-        '--minStart', type=int, default=None,
+        '--minStart', type=int,
         help='Reads that start before this subject offset should not be '
         'shown.')
 
     parser.add_argument(
-        '--maxStop', type=int, default=None,
+        '--maxStop', type=int,
         help='Reads that end after this subject offset should not be shown.')
 
     parser.add_argument(
@@ -128,94 +128,99 @@ if __name__ == '__main__':
         help='If True, only keep the best alignment for each read.')
 
     parser.add_argument(
-        '--maxAlignmentsPerRead', type=int, default=None,
+        '--maxAlignmentsPerRead', type=int,
         help=('Reads with more than this many alignments will be elided. Pass '
               'zero to only keep reads with no matches (alignments).'))
 
     parser.add_argument(
-        '--scoreCutoff', type=float, default=None,
+        '--scoreCutoff', type=float,
         help=('A float score. Matches with scores worse than this will be '
               'ignored.'))
 
     parser.add_argument(
-        '--maxHspsPerHit', type=int, default=None,
+        '--percentageIdenticalCutoff', type=float,
+        help=('A float percentage identity (0.0 to 100.0). Matches with '
+              'percent identity scores less than this will be ignored.'))
+
+    parser.add_argument(
+        '--maxHspsPerHit', type=int,
         help='A numeric max number of HSPs to show for each hit on hitId.')
 
     parser.add_argument(
-        '--whitelist', nargs='+', default=None, action='append',
+        '--whitelist', nargs='+', action='append',
         help='sequence titles that should be whitelisted')
 
     parser.add_argument(
-        '--whitelistFile', default=None,
+        '--whitelistFile',
         help=('A file containing sequence titles (one per line) to be '
               'whitelisted'))
 
     parser.add_argument(
-        '--blacklist', nargs='+', default=None, action='append',
+        '--blacklist', nargs='+', action='append',
         help='sequence titles that should be blacklisted')
 
     parser.add_argument(
-        '--blacklistFile', default=None,
+        '--blacklistFile',
         help=('A file containing sequence titles (one per line) to be '
               'blacklisted'))
 
     parser.add_argument(
-        '--titleRegex', default=None,
+        '--titleRegex',
         help='a regex that sequence titles must match.')
 
     parser.add_argument(
-        '--negativeTitleRegex', default=None,
+        '--negativeTitleRegex',
         help='a regex that sequence titles must not match.')
 
     parser.add_argument(
-        '--truncateTitlesAfter', default=None,
+        '--truncateTitlesAfter',
         help=('a string that titles will be truncated beyond. If the '
               'truncated version of a title has already been seen, '
               'that title will be skipped.'))
 
     parser.add_argument(
-        '--minSequenceLen', type=int, default=None,
+        '--minSequenceLen', type=int,
         help='sequences of lesser length will be elided.')
 
     parser.add_argument(
-        '--maxSequenceLen', type=int, default=None,
+        '--maxSequenceLen', type=int,
         help='sequences of greater length will be elided.')
 
     parser.add_argument(
-        '--taxonomy', default=None,
+        '--taxonomy',
         help=('a string of the taxonomic group on which should be '
               'filtered. eg "Vira" will filter on viruses.'))
 
     # Args for filtering on TitlesAlignments.
     parser.add_argument(
-        '--minMatchingReads', type=int, default=None,
+        '--minMatchingReads', type=int,
         help='sequences that are matched by fewer reads will be elided.')
 
     parser.add_argument(
-        '--minMedianScore', type=float, default=None,
+        '--minMedianScore', type=float,
         help=('sequences that are matched with a median score that is '
               'worse will be elided.'))
 
     parser.add_argument(
-        '--withScoreBetterThan', type=float, default=None,
+        '--withScoreBetterThan', type=float,
         help=('sequences that are matched without at least one score '
               'at least this good will be elided.'))
 
     parser.add_argument(
-        '--minNewReads', type=float, default=None,
+        '--minNewReads', type=float,
         help=('The fraction of its reads by which a new read set must differ '
               'from all previously seen read sets in order to be considered '
               'acceptably different.'))
 
     parser.add_argument(
-        '--maxTitles', type=int, default=None,
+        '--maxTitles', type=int,
         help=('The maximum number of titles to keep. If more titles than '
               'this result from the filtering, titles will be sorted '
               '(according to the --sortOn value) and only the best will be '
               'retained.'))
 
     parser.add_argument(
-        '--minCoverage', type=float, default=None,
+        '--minCoverage', type=float,
         help=('The (0.0 to 1.0) minimum fraction of a subject sequence that '
               'must be matched by at least one read.'))
 
@@ -231,7 +236,7 @@ if __name__ == '__main__':
               'the score.'))
 
     parser.add_argument(
-        '--outputDir', default=None, required=True,
+        '--outputDir', required=True,
         help='Specifies a directory to write the HTML summary to. Required.')
 
     parser.add_argument(
@@ -334,6 +339,7 @@ if __name__ == '__main__':
         oneAlignmentPerRead=args.oneAlignmentPerRead,
         maxHspsPerHit=args.maxHspsPerHit,
         scoreCutoff=args.scoreCutoff,
+        percentageIdenticalCutoff=args.percentageIdenticalCutoff,
         whitelist=whitelist, blacklist=blacklist,
         whitelistFile=args.whitelistFile, blacklistFile=args.blacklistFile,
         titleRegex=args.titleRegex, negativeTitleRegex=args.negativeTitleRegex,
diff --git a/bin/noninteractive-alignment-panel.py b/bin/noninteractive-alignment-panel.py
index 44152fc..77a4297 100755
--- a/bin/noninteractive-alignment-panel.py
+++ b/bin/noninteractive-alignment-panel.py
@@ -134,12 +134,12 @@ if __name__ == '__main__':
 
     # Args for filtering on ReadsAlignments.
     parser.add_argument(
-        '--minStart', type=int, default=None,
+        '--minStart', type=int,
         help='Reads that start before this subject offset should not be '
         'shown.')
 
     parser.add_argument(
-        '--maxStop', type=int, default=None,
+        '--maxStop', type=int,
         help='Reads that end after this subject offset should not be shown.')
 
     parser.add_argument(
@@ -147,94 +147,99 @@ if __name__ == '__main__':
         help='If True, only keep the best alignment for each read.')
 
     parser.add_argument(
-        '--maxAlignmentsPerRead', type=int, default=None,
+        '--maxAlignmentsPerRead', type=int,
         help=('Reads with more than this many alignments will be elided. Pass '
               'zero to only keep reads with no matches (alignments).'))
 
     parser.add_argument(
-        '--scoreCutoff', type=float, default=None,
+        '--scoreCutoff', type=float,
         help=('A float score. Matches with scores worse than this will be '
               'ignored.'))
 
     parser.add_argument(
-        '--maxHspsPerHit', type=int, default=None,
+        '--percentageIdenticalCutoff', type=float,
+        help=('A float percentage identity (0.0 to 100.0). Matches with '
+              'percent identity scores less than this will be ignored.'))
+
+    parser.add_argument(
+        '--maxHspsPerHit', type=int,
         help='A numeric max number of HSPs to show for each hit on hitId.')
 
     parser.add_argument(
-        '--whitelist', nargs='+', default=None, action='append',
+        '--whitelist', nargs='+', action='append',
         help='sequence titles that should be whitelisted')
 
     parser.add_argument(
-        '--whitelistFile', default=None,
+        '--whitelistFile',
         help=('A file containing sequence titles (one per line) to be '
               'whitelisted'))
 
     parser.add_argument(
-        '--blacklist', nargs='+', default=None, action='append',
+        '--blacklist', nargs='+', action='append',
         help='sequence titles that should be blacklisted')
 
     parser.add_argument(
-        '--blacklistFile', default=None,
+        '--blacklistFile',
         help=('A file containing sequence titles (one per line) to be '
               'blacklisted'))
 
     parser.add_argument(
-        '--titleRegex', default=None,
+        '--titleRegex',
         help='a regex that sequence titles must match.')
 
     parser.add_argument(
-        '--negativeTitleRegex', default=None,
+        '--negativeTitleRegex',
         help='a regex that sequence titles must not match.')
 
     parser.add_argument(
-        '--truncateTitlesAfter', default=None,
+        '--truncateTitlesAfter',
         help=('a string that titles will be truncated beyond. If the '
               'truncated version of a title has already been seen, '
               'that title will be skipped.'))
 
     parser.add_argument(
-        '--minSequenceLen', type=int, default=None,
+        '--minSequenceLen', type=int,
         help='sequences of lesser length will be elided.')
 
     parser.add_argument(
-        '--maxSequenceLen', type=int, default=None,
+        '--maxSequenceLen', type=int,
         help='sequences of greater length will be elided.')
 
     parser.add_argument(
-        '--taxonomy', default=None,
+        '--taxonomy',
         help=('a string of the taxonomic group on which should be '
               'filtered. eg "Vira" will filter on viruses.'))
 
     # Args for filtering on TitlesAlignments.
     parser.add_argument(
-        '--minMatchingReads', type=int, default=None,
+        '--minMatchingReads', type=int,
         help='sequences that are matched by fewer reads will be elided.')
 
     parser.add_argument(
-        '--minMedianScore', type=float, default=None,
+        '--minMedianScore', type=float,
         help=('sequences that are matched with a median score that is '
               'worse will be elided.'))
 
     parser.add_argument(
-        '--withScoreBetterThan', type=float, default=None,
+        '--withScoreBetterThan', type=float,
         help=('sequences that are matched without at least one score '
               'at least this good will be elided.'))
 
     parser.add_argument(
-        '--minNewReads', type=float, default=None,
+        '--minNewReads', type=float,
         help=('The fraction of its reads by which a new read set must differ '
               'from all previously seen read sets in order to be considered '
               'acceptably different.'))
 
     parser.add_argument(
-        '--maxTitles', type=int, default=None,
+        '--maxTitles', type=int,
         help=('The maximum number of titles to keep. If more titles than '
               'this result from the filtering, titles will be sorted '
               '(according to the --sortOn value) and only the best will be '
               'retained.'))
 
     parser.add_argument(
-        '--minCoverage', type=float, default=None,
+        '--minCoverage', type=float,
         help=('The (0.0 to 1.0) minimum fraction of a subject sequence that '
               'must be matched by at least one read.'))
 
@@ -250,7 +255,7 @@ if __name__ == '__main__':
               'the score.'))
 
     parser.add_argument(
-        '--outputDir', default=None, required=True,
+        '--outputDir', required=True,
         help='Specifies a directory to write the HTML summary to. Required.')
 
     parser.add_argument(
@@ -384,6 +389,7 @@ if __name__ == '__main__':
         oneAlignmentPerRead=args.oneAlignmentPerRead,
         maxHspsPerHit=args.maxHspsPerHit,
         scoreCutoff=args.scoreCutoff,
+        percentageIdenticalCutoff=args.percentageIdenticalCutoff,
         whitelist=whitelist, blacklist=blacklist,
         whitelistFile=args.whitelistFile, blacklistFile=args.blacklistFile,
         titleRegex=args.titleRegex, negativeTitleRegex=args.negativeTitleRegex,
diff --git a/dark/alignments.py b/dark/alignments.py
index 6ada96c..89aa032 100644
--- a/dark/alignments.py
+++ b/dark/alignments.py
@@ -105,6 +105,9 @@ class ReadsAlignmentsFilter(object):
         alignment for each read.
     @param scoreCutoff: A C{float} score. Matches with scores that are not
         better than this score will be ignored.
+    @param percentageIdenticalCutoff: A C{float} percentage identity (i.e.,
+        0.0 to 100.0 NOT a fraction). Matches with percentage identity less
+        than this will be ignored.
     @param whitelist: If not C{None}, a set of exact titles that are always
         acceptable (though the match info for a whitelist title may rule it
         out for other reasons).
@@ -130,7 +133,8 @@ class ReadsAlignmentsFilter(object):
                  minSequenceLen=None, maxSequenceLen=None,
                  minStart=None, maxStop=None,
                  oneAlignmentPerRead=False, maxHspsPerHit=None,
-                 scoreCutoff=None, whitelist=None, blacklist=None,
+                 scoreCutoff=None, percentageIdenticalCutoff=None,
+                 whitelist=None, blacklist=None,
                  whitelistFile=None, blacklistFile=None,
                  titleRegex=None, negativeTitleRegex=None,
                  truncateTitlesAfter=None, taxonomy=None, readIdRegex=None):
@@ -144,6 +148,7 @@ class ReadsAlignmentsFilter(object):
         self.oneAlignmentPerRead = oneAlignmentPerRead
         self.maxHspsPerHit = maxHspsPerHit
         self.scoreCutoff = scoreCutoff
+        self.percentageIdenticalCutoff = percentageIdenticalCutoff
 
         # If we've been asked to filter on matched sequence titles in any way,
         # build a title filter.
@@ -297,6 +302,29 @@ class ReadsAlignmentsFilter(object):
             else:
                 return False
 
+        # Throw out HSPs whose percentage identical is not good enough.
+        #
+        # Note that if the percentIdentical is None it means that we did
+        # not run DIAMOND or BLAST in a way that generated that output, so
+        # we cannot eliminate such an HSP.
+        if self.percentageIdenticalCutoff is not None:
+            cutoff = self.percentageIdenticalCutoff
+            wantedAlignments = []
+            for alignment in readAlignments:
+                hsps = alignment.hsps
+                wantedHsps = []
+                for hsp in hsps:
+                    percentIdentical = hsp.percentIdentical
+                    if percentIdentical is None or percentIdentical >= cutoff:
+                        wantedHsps.append(hsp)
+                if wantedHsps:
+                    alignment.hsps = wantedHsps
+                    wantedAlignments.append(alignment)
+            if wantedAlignments:
+                readAlignments[:] = wantedAlignments
+            else:
+                return False
+
         # Throw out HSPs that don't match in the desired place on the
         # matched sequence.
         minStart = self.minStart
diff --git a/dark/civ/proteins.py b/dark/civ/proteins.py
index f81d546..fa6e6b0 100644
--- a/dark/civ/proteins.py
+++ b/dark/civ/proteins.py
@@ -613,7 +613,6 @@ class ProteinGrouper(object):
                 '<li>All read lengths (in parentheses).</li>')
 
         proteinFieldsDescription.extend([
-            '<li>Protein name.</li>',
             '</ol>',
             '</p>',
         ])
@@ -669,7 +668,8 @@ class ProteinGrouper(object):
         result.extend(proteinFieldsDescription)
 
         # Write a linked table of contents by pathogen.
-        append('<p><span class="index-name">Pathogen index:</span>')
+        append('<p><span class="index-name">%s index:</span>' % (
+            'Bacterium' if pathogenType == 'bacterial' else 'Virus'))
         append('<span class="index">')
         for genomeAccession in genomeAccessions:
             genomeInfo = self._db.findGenome(genomeAccession)
@@ -807,7 +807,7 @@ class ProteinGrouper(object):
                         countClass = readCountColors.thresholdToCssName(
                             readCountColors.thresholdForCount(
                                 proteinMatch['readCount']))
-                        appendNoSpace('<span class="%s">%3s</span>' % (
+                        appendNoSpace('<span class="%s">%4s</span>' % (
                             countClass, proteinMatch['readAndHspCountStr']))
                     else:
                         appendNoSpace('%(readAndHspCountStr)3s' % proteinMatch)
@@ -861,13 +861,15 @@ class ProteinGrouper(object):
         if bootstrapTreeviewDir:
             append('''
                 <script>
-                  var tree = %s;
-                  $('#tree').treeview({
-                      data: tree,
-                      enableLinks: true,
-                      levels: 0,
-                  });
-               </script>
+                $(document).ready(function(){
+                    var tree = %s;
+                    $('#tree').treeview({
+                        data: tree,
+                        enableLinks: true,
+                        levels: 0,
+                    });
+                });
+                </script>
             ''' % taxonomyHierarchy.toJSON())
 
         # Write all samples (with pathogens (with proteins)).
@@ -945,7 +947,7 @@ class ProteinGrouper(object):
                         countClass = readCountColors.thresholdToCssName(
                             readCountColors.thresholdForCount(
                                 proteinMatch['readCount']))
-                        appendNoSpace('<span class="%s">%3s</span>' % (
+                        appendNoSpace('<span class="%s">%4s</span>' % (
                             countClass, proteinMatch['readAndHspCountStr']))
                     else:
                         appendNoSpace('%(readAndHspCountStr)3s' % proteinMatch)
@@ -1025,14 +1027,17 @@ class ProteinGrouper(object):
         ax.set_xticks([])
         ax.set_xlim((-0.5, nSamples - 0.5))
         ax.vlines(x, yMin, readCounts, color=color)
+
+        genomeInfo = self._db.findGenome(genomeAccession)
+
         if highlighted:
             title = '%s\nIn red: %s' % (
-                genomeAccession, fill(', '.join(highlighted), 50))
+                genomeInfo['organism'], fill(', '.join(highlighted), 50))
         else:
             # Add a newline to keep the first line of each title at the
             # same place as those titles that have an "In red:" second
             # line.
-            title = genomeAccession + '\n'
+            title = genomeInfo['organism'] + '\n'
 
         ax.set_title(title, fontsize=10)
         ax.tick_params(axis='both', which='major', labelsize=8)
diff --git a/dark/diamond/conversion.py b/dark/diamond/conversion.py
index 234545a..4e0ecdf 100644
--- a/dark/diamond/conversion.py
+++ b/dark/diamond/conversion.py
@@ -123,12 +123,25 @@ class DiamondTabularFormatReader(object):
     Provide a method that yields parsed tabular records from a file. Store and
     make accessible the global DIAMOND parameters.
 
-    Make sure you run DIAMOND with the right output format. You must use:
-        --outfmt 6 qtitle stitle bitscore evalue qframe qseq qstart qend sseq
-                   sstart send slen btop nident positive
+    Make sure you run DIAMOND with the right output format. You must have run
+    with --outfmt 6 and one of the following sets of arguments (in the given
+    order):
+
+        13 arg version:
+        qtitle stitle bitscore evalue qframe qseq qstart qend sseq sstart send
+            slen btop
+
+        14 arg version:
+        qtitle stitle bitscore evalue qframe qseq qstart qend sseq sstart send
+            slen btop pident
+
+        15 arg version:
+        qtitle stitle bitscore evalue qframe qseq qstart qend sseq sstart send
+            slen btop nident positive
 
     @param filename: A C{str} filename or an open file pointer, containing
         DIAMOND tabular records.
+
     """
 
     def __init__(self, filename):
@@ -159,36 +172,36 @@ class DiamondTabularFormatReader(object):
             record = {}
             for line in fp:
                 line = line[:-1]
-                try:
-                    (qtitle, stitle, bitscore, evalue, qframe, qseq,
-                     qstart, qend, sseq, sstart, send, slen, btop, nident,
-                     positive) = line.split('\t')
-                except ValueError as e:
-                    # We may not be able to find 'nident' and 'positives'
-                    # because they were added in version 2.0.3 and will not
-                    # be present in any of our JSON output generated before
-                    # that. So those values will be None when reading
-                    # DIAMOND output without those fields, but that's much
-                    # better than no longer being able to read that data.
-                    if six.PY2:
-                        error = 'need more than 13 values to unpack'
-                    else:
-                        error = (
-                            'not enough values to unpack (expected 15, '
-                            'got 13)')
-                    if str(e) == error:
-                        (qtitle, stitle, bitscore, evalue, qframe,
-                         qseq, qstart, qend, sseq, sstart, send, slen,
-                         btop) = line.split('\t')
-                        nident = positive = None
-                    else:
-                        raise
+                fields = line.split('\t')
+                nFields = len(fields)
+
+                if nFields == 13:
+                    (qtitle, stitle, bitscore, evalue, qframe, qseq, qstart,
+                     qend, sseq, sstart, send, slen, btop) = fields
+                    nident = positive = pident = None
+                elif nFields == 14:
+                    (qtitle, stitle, bitscore, evalue, qframe, qseq, qstart,
+                     qend, sseq, sstart, send, slen, btop, pident) = fields
+                    nident = positive = None
+                elif nFields == 15:
+                    (qtitle, stitle, bitscore, evalue, qframe, qseq, qstart,
+                     qend, sseq, sstart, send, slen, btop, nident,
+                     positive) = fields
+                    pident = None
+                else:
+                    raise ValueError(
+                        'Could not make sense of DIAMOND output. You must use '
+                        '--outfmt 6 with 13, 14, or 15 arguments. See %s for '
+                        'their names.' % __file__)
+
                 hsp = {
                     'bits': float(bitscore),
                     'btop': btop,
                     'expect': float(evalue),
                     'frame': int(qframe),
                     'identicalCount': None if nident is None else int(nident),
+                    'percentIdentical': (
+                        None if pident is None else float(pident)),
                     'positiveCount': (
                         None if positive is None else int(positive)),
                     'query': qseq,
@@ -198,6 +211,7 @@ class DiamondTabularFormatReader(object):
                     'sbjct_start': int(sstart),
                     'sbjct_end': int(send),
                 }
+
                 if previousQtitle == qtitle:
                     # We have already started accumulating alignments for this
                     # query.
@@ -346,14 +360,17 @@ class JSONRecordsReader(object):
                     subjectEnd=normalized['subjectEnd'],
                     readMatchedSequence=diamondHsp['query'],
                     subjectMatchedSequence=diamondHsp['sbjct'],
-                    # Use blastHsp.get on identicalCount and positiveCount
-                    # because they were added in version 2.0.3 and will not
-                    # be present in any of our JSON output generated before
-                    # that. Those values will be None for those JSON files,
-                    # but that's much better than no longer being able to
-                    # read all that data.
+                    # Use diamondHsp.get on identicalCount, positiveCount,
+                    # and percentIdentical because they were either added
+                    # in version 2.0.3 or we didn't start using them until
+                    # much later and so will not be present in any of our
+                    # JSON output generated before that. Those values will
+                    # be None when reading those JSON files, but that's
+                    # much better than no longer being able to read all
+                    # that earlier data.
                     identicalCount=diamondHsp.get('identicalCount'),
-                    positiveCount=diamondHsp.get('positiveCount'))
+                    positiveCount=diamondHsp.get('positiveCount'),
+                    percentIdentical=diamondHsp.get('percentIdentical'))
 
                 alignment.addHsp(hsp)
 
diff --git a/dark/hsp.py b/dark/hsp.py
index f25d381..742108c 100644
--- a/dark/hsp.py
+++ b/dark/hsp.py
@@ -62,12 +62,14 @@ class _Base(object):
         and query had a positive score in the scoring matrix used during
         matching (this is probably only different from the C{identicalCount}
         when matching amino acids (i.e., not nucleotides).
+    @param percentIdentical: A C{float} percentage (i.e., ranging from 0.0 to
+        100.0, NOT a fraction) of amino acids that were identical in the match.
     """
     def __init__(self, readStart=None, readEnd=None, readStartInSubject=None,
                  readEndInSubject=None, readFrame=None, subjectStart=None,
                  subjectEnd=None, subjectFrame=None, readMatchedSequence=None,
                  subjectMatchedSequence=None, identicalCount=None,
-                 positiveCount=None):
+                 positiveCount=None, percentIdentical=None):
         self.readStart = readStart
         self.readEnd = readEnd
         self.readStartInSubject = readStartInSubject
@@ -80,6 +82,7 @@ class _Base(object):
         self.subjectMatchedSequence = subjectMatchedSequence
         self.identicalCount = identicalCount
         self.positiveCount = positiveCount
+        self.percentIdentical = percentIdentical
 
     def __lt__(self, other):
         return self.score < other.score
@@ -115,6 +118,7 @@ class _Base(object):
             'subjectMatchedSequence': self.subjectMatchedSequence,
             'identicalCount': self.identicalCount,
             'positiveCount': self.positiveCount,
+            'percentIdentical': self.percentIdentical,
         }
 
 


Test patch:
diff --git a/test/blast/test_titles.py b/test/blast/test_titles.py
index 70cded1..b4aa345 100644
--- a/test/blast/test_titles.py
+++ b/test/blast/test_titles.py
@@ -313,6 +313,7 @@ class TestTitlesAlignments(TestCase):
                                 'hsps': [
                                     {
                                         'identicalCount': None,
+                                        'percentIdentical': None,
                                         'positiveCount': None,
                                         'readEnd': 68,
                                         'readEndInSubject': 1405,
@@ -337,6 +338,7 @@ class TestTitlesAlignments(TestCase):
                                 'hsps': [
                                     {
                                         'identicalCount': None,
+                                        'percentIdentical': None,
                                         'positiveCount': None,
                                         'readEnd': 68,
                                         'readEndInSubject': 1405,
@@ -368,6 +370,7 @@ class TestTitlesAlignments(TestCase):
                                 'hsps': [
                                     {
                                         'identicalCount': None,
+                                        'percentIdentical': None,
                                         'positiveCount': None,
                                         'readEnd': 68,
                                         'readEndInSubject': 11405,
@@ -399,6 +402,7 @@ class TestTitlesAlignments(TestCase):
                                 'hsps': [
                                     {
                                         'identicalCount': None,
+                                        'percentIdentical': None,
                                         'positiveCount': None,
                                         'readEnd': 68,
                                         'readEndInSubject': 10405,
@@ -430,6 +434,7 @@ class TestTitlesAlignments(TestCase):
                                 'hsps': [
                                     {
                                         'identicalCount': None,
+                                        'percentIdentical': None,
                                         'positiveCount': None,
                                         'readEnd': 68,
                                         'readEndInSubject': 15405,
@@ -461,6 +466,7 @@ class TestTitlesAlignments(TestCase):
                                 'hsps': [
                                     {
                                         'identicalCount': None,
+                                        'percentIdentical': None,
                                         'positiveCount': None,
                                         'readEnd': 68,
                                         'readEndInSubject': 12405,
diff --git a/test/diamond/test_alignments.py b/test/diamond/test_alignments.py
index 3245c6e..9ac2e4e 100644
--- a/test/diamond/test_alignments.py
+++ b/test/diamond/test_alignments.py
@@ -783,6 +783,73 @@ class TestDiamondReadsAlignmentsFiltering(TestCase):
             self.assertEqual(1, len(result[0]))
             self.assertEqual('Merkel2', result[0][0].subjectTitle)
 
+    def testPercentIdentityCutoffRemovesEntireAlignment(self):
+        """
+        If the L{DiamondReadsAlignments} filter function is supposed to filter
+        on a percent identity cut-off and the cut-off value results in an
+        alignment with no HSPs, then the alignment must be removed entirely.
+        """
+        record = {
+            "query": "H6E8I1T01BFUH9",
+            "alignments": [
+                {
+                    "length": 961,
+                    "hsps": [
+                        {
+                            "sbjct_end": 869,
+                            "expect": 1.25854e-10,
+                            "sbjct": "AAAAAAAAAAAAAAAAAAAAAAAAAAAAAA",
+                            "sbjct_start": 836,
+                            "query": ("AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA"
+                                      "AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA"
+                                      "AAAAAAAAAAAAAAAAAAAAAAAAAA"),
+                            "frame": 1,
+                            "query_end": 462,
+                            "bits": 150,
+                            "btop": "",
+                            "query_start": 362,
+                            "percentIdentical": 40.0,
+                        },
+                    ],
+                    "title": "Merkel1"
+                },
+                {
+                    "length": 740,
+                    "hsps": [
+                        {
+                            "sbjct_end": 647,
+                            "expect": 1.25e-43,
+                            "sbjct": "AAAAAAAAAAAAAAAAAAAAAAAAAAAAAA",
+                            "sbjct_start": 614,
+                            "query": ("AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA"
+                                      "AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA"
+                                      "AAAAAAAAAAAAAAAAAAAAAAAAAA"),
+                            "frame": 1,
+                            "query_end": 462,
+                            "bits": 180,
+                            "btop": "",
+                            "query_start": 362,
+                            "percentIdentical": 45.0,
+                        }
+                    ],
+                    "title": "Merkel2"
+                }
+            ]
+        }
+
+        mockOpener = mockOpen(read_data=dumps(PARAMS) + '\n' +
+                              dumps(record) + '\n')
+        with patch.object(builtins, 'open', mockOpener):
+            reads = Reads()
+            reads.add(Read('H6E8I1T01BFUH9', 'A' * 500))
+            readsAlignments = DiamondReadsAlignments(
+                reads, 'file.json', databaseFilename='database.fasta')
+            result = list(readsAlignments.filter(
+                percentageIdenticalCutoff=41.0))
+            self.assertEqual(1, len(result))
+            self.assertEqual(1, len(result[0]))
+            self.assertEqual('Merkel2', result[0][0].subjectTitle)
+
     def testScoreCutoffRemovesHsps_Bits(self):
         """
         If the L{DiamondRecords} records function is supposed to filter on
@@ -952,6 +1019,182 @@ class TestDiamondReadsAlignmentsFiltering(TestCase):
             self.assertEqual(1, len(result[0][1].hsps))
             self.assertEqual(LSP(1.25e-43), result[0][1].hsps[0])
 
+    def testPercentIdentityCutoffRemovesHsps(self):
+        """
+        If the L{DiamondRecords} records function is supposed to filter on
+        percentage identical cut-off and the cut-off value results in some HSPs
+        being invalid, then those HSPs must be removed entirely.
+        """
+        record = {
+            "query": "H6E8I1T01BFUH9",
+            "alignments": [
+                {
+                    "length": 961,
+                    "hsps": [
+                        {
+                            "sbjct_end": 869,
+                            "expect": 1.25854e-10,
+                            "sbjct": "AAAAAAAAAAAAAAAAAAAAAAAAAAAAAA",
+                            "sbjct_start": 836,
+                            "query": ("AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA"
+                                      "AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA"
+                                      "AAAAAAAAAAAAAAAAAAAAAAAAAA"),
+                            "frame": 1,
+                            "query_end": 462,
+                            "bits": 150,
+                            "btop": "",
+                            "query_start": 362,
+                            "percentIdentical": 30.0,
+                        },
+                        {
+                            "sbjct_end": 869,
+                            "expect": 1.25e-20,
+                            "sbjct": "AAAAAAAAAAAAAAAAAAAAAAAAAAAAAA",
+                            "sbjct_start": 836,
+                            "query": ("AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA"
+                                      "AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA"
+                                      "AAAAAAAAAAAAAAAAAAAAAAAAAA"),
+                            "frame": 1,
+                            "query_end": 462,
+                            "bits": 170,
+                            "btop": "",
+                            "query_start": 362,
+                            "percentIdentical": 70.0,
+                        }
+                    ],
+                    "title": "Merkel1"
+                },
+                {
+                    "length": 740,
+                    "hsps": [
+                        {
+                            "sbjct_end": 647,
+                            "expect": 1.25e-43,
+                            "sbjct": "AAAAAAAAAAAAAAAAAAAAAAAAAAAAAA",
+                            "sbjct_start": 614,
+                            "query": ("AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA"
+                                      "AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA"
+                                      "AAAAAAAAAAAAAAAAAAAAAAAAAA"),
+                            "frame": 1,
+                            "query_end": 462,
+                            "bits": 180,
+                            "btop": "",
+                            "query_start": 362,
+                            "percentIdentical": 80.0,
+                        }
+                    ],
+                    "title": "Merkel2"
+                }
+            ]
+        }
+
+        mockOpener = mockOpen(read_data=dumps(PARAMS) + '\n' +
+                              dumps(record) + '\n')
+        with patch.object(builtins, 'open', mockOpener):
+            reads = Reads()
+            reads.add(Read('H6E8I1T01BFUH9', 'A' * 500))
+            readsAlignments = DiamondReadsAlignments(
+                reads, 'file.json', databaseFilename='database.fasta')
+            result = list(readsAlignments.filter(
+                percentageIdenticalCutoff=50.0))
+
+            # There should only be one HSP left in the alignments for the
+            # first read, and it should have the right score (bit score).
+            self.assertEqual(1, len(result[0][0].hsps))
+            self.assertEqual(70.0, result[0][0].hsps[0].percentIdentical)
+
+            # The second alignment should also be present.
+            self.assertEqual(1, len(result[0][1].hsps))
+            self.assertEqual(80.0, result[0][1].hsps[0].percentIdentical)
+
+    def testPercentIdentityCutoffIgnoresNoneValues(self):
+        """
+        If the L{DiamondRecords} records function is supposed to filter on
+        percentage identical cut-off and the cut-off value in some HSPs
+        are None, then those HSPs must not be removed.
+        """
+        record = {
+            "query": "H6E8I1T01BFUH9",
+            "alignments": [
+                {
+                    "length": 961,
+                    "hsps": [
+                        {
+                            "sbjct_end": 869,
+                            "expect": 1.25854e-10,
+                            "sbjct": "AAAAAAAAAAAAAAAAAAAAAAAAAAAAAA",
+                            "sbjct_start": 836,
+                            "query": ("AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA"
+                                      "AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA"
+                                      "AAAAAAAAAAAAAAAAAAAAAAAAAA"),
+                            "frame": 1,
+                            "query_end": 462,
+                            "bits": 150,
+                            "btop": "",
+                            "query_start": 362,
+                            "percentIdentical": None,
+                        },
+                        {
+                            "sbjct_end": 869,
+                            "expect": 1.25e-20,
+                            "sbjct": "AAAAAAAAAAAAAAAAAAAAAAAAAAAAAA",
+                            "sbjct_start": 836,
+                            "query": ("AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA"
+                                      "AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA"
+                                      "AAAAAAAAAAAAAAAAAAAAAAAAAA"),
+                            "frame": 1,
+                            "query_end": 462,
+                            "bits": 170,
+                            "btop": "",
+                            "query_start": 362,
+                            "percentIdentical": 70.0,
+                        }
+                    ],
+                    "title": "Merkel1"
+                },
+                {
+                    "length": 740,
+                    "hsps": [
+                        {
+                            "sbjct_end": 647,
+                            "expect": 1.25e-43,
+                            "sbjct": "AAAAAAAAAAAAAAAAAAAAAAAAAAAAAA",
+                            "sbjct_start": 614,
+                            "query": ("AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA"
+                                      "AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA"
+                                      "AAAAAAAAAAAAAAAAAAAAAAAAAA"),
+                            "frame": 1,
+                            "query_end": 462,
+                            "bits": 180,
+                            "btop": "",
+                            "query_start": 362,
+                            "percentIdentical": 80.0,
+                        }
+                    ],
+                    "title": "Merkel2"
+                }
+            ]
+        }
+
+        mockOpener = mockOpen(read_data=dumps(PARAMS) + '\n' +
+                              dumps(record) + '\n')
+        with patch.object(builtins, 'open', mockOpener):
+            reads = Reads()
+            reads.add(Read('H6E8I1T01BFUH9', 'A' * 500))
+            readsAlignments = DiamondReadsAlignments(
+                reads, 'file.json', databaseFilename='database.fasta')
+            result = list(readsAlignments.filter(
+                percentageIdenticalCutoff=75.0))
+
+            # Only one HSP should be left in the alignments for the first
+            # read, and it should have a None percent identity score.
+            self.assertEqual(1, len(result[0][0].hsps))
+            self.assertIs(None, result[0][0].hsps[0].percentIdentical)
+
+            # The second alignment should also be present.
+            self.assertEqual(1, len(result[0][1].hsps))
+            self.assertEqual(80.0, result[0][1].hsps[0].percentIdentical)
+
     def testTitleByRegexCaseInvariant(self):
         """
         Filtering with a title regex must work independent of case.
diff --git a/test/diamond/test_conversion.py b/test/diamond/test_conversion.py
index 93cad66..4238b35 100644
--- a/test/diamond/test_conversion.py
+++ b/test/diamond/test_conversion.py
@@ -64,6 +64,28 @@ BHAV	TAIV	28.1	0.008	1	PKELHGLI	14	118	SLKSKE	15	131	307	8
 BHAV	SouthBay	28.1	0.009	1	CRPTF	4	293	EFVFIY	6	342	343	5
 """
 
+# The 14 fields expected in the DIAMOND output with pident (percent identical).
+#
+# qtitle, stitle, bitscore, evalue, qframe, qseq, qstart, qend, sseq, sstart,
+# send, slen, btop, pident
+#
+# See the --outfmt section of 'diamond help' for detail on these directives.
+#
+# Note that the fields below must be separated by TABs.
+DIAMOND_RECORDS_WITH_PERCENT_IDENTICAL = """\
+ACC94	INSV	29.6	0.003	1	EFII	178	295	SSSEV	175	285	295	4	5.0
+ACC94	CASV	28.1	0.008	1	KLL	7	37	ITRV	9	39	300	3	10.0
+ACC94	GoldenGate	28.1	0.009	1	IKSKL	7	35	EETSR	9	37	293	5	15.0
+ACC94	GoldenGate	23.5	0.21	1	TIMSVV	177	240	DDMV	179	235	293	6	20.0
+ACC94	InfluenzaC	25.0	0.084	1	LHVNYL	1	203	DEELKA	2	210	290	6	25.0
+ACC94	InfluenzaC	18.5	9.1	1	SEIICEVLK	226	257	VETVAQ	20	45	290	9	30.0
+ACC94	FERV	24.6	0.11	1	YSCFT-NSEK	176	276	LGKRMFC	152	243	270	10	35.0
+AKAV	AKAV	634	0.0	1	GEPFSVYG	1	306	NIYGEP	1	306	306	8	40.0
+AKAV	WYOV	401	7e-143	1	PFSVYGRF	1	306	GEPMS	1	294	294	8	45.0
+BHAV	TAIV	28.1	0.008	1	PKELHGLI	14	118	SLKSKE	15	131	307	8	50.0
+BHAV	SouthBay	28.1	0.009	1	CRPTF	4	293	EFVFIY	6	342	343	5	55.0
+"""
+
 DIAMOND_RECORD_WITH_SPACES_IN_TITLES = """\
 ACC 94	IN SV	29.6	0.003	1	EFII	178	295	SSSEV	175	285	295	4
 """
@@ -87,6 +109,7 @@ DIAMOND_RECORDS_DUMPED = '\n'.join([
                         "expect": 0.003,
                         "frame": 1,
                         "identicalCount": 0,
+                        "percentIdentical": None,
                         "positiveCount": 1,
                         "query": "EFII",
                         "query_end": 295,
@@ -107,6 +130,7 @@ DIAMOND_RECORDS_DUMPED = '\n'.join([
                         "expect": 0.008,
                         "frame": 1,
                         "identicalCount": 1,
+                        "percentIdentical": None,
                         "positiveCount": 2,
                         "query": "KLL",
                         "query_end": 37,
@@ -127,6 +151,7 @@ DIAMOND_RECORDS_DUMPED = '\n'.join([
                         "expect": 0.009,
                         "frame": 1,
                         "identicalCount": 2,
+                        "percentIdentical": None,
                         "positiveCount": 3,
                         "query": "IKSKL",
                         "query_end": 35,
@@ -141,6 +166,7 @@ DIAMOND_RECORDS_DUMPED = '\n'.join([
                         "expect": 0.21,
                         "frame": 1,
                         "identicalCount": 3,
+                        "percentIdentical": None,
                         "positiveCount": 4,
                         "query": "TIMSVV",
                         "query_end": 240,
@@ -161,6 +187,7 @@ DIAMOND_RECORDS_DUMPED = '\n'.join([
                         "expect": 0.084,
                         "frame": 1,
                         "identicalCount": 4,
+                        "percentIdentical": None,
                         "positiveCount": 5,
                         "query": "LHVNYL",
                         "query_end": 203,
@@ -175,6 +202,7 @@ DIAMOND_RECORDS_DUMPED = '\n'.join([
                         "expect": 9.1,
                         "frame": 1,
                         "identicalCount": 5,
+                        "percentIdentical": None,
                         "positiveCount": 6,
                         "query": "SEIICEVLK",
                         "query_end": 257,
@@ -195,6 +223,7 @@ DIAMOND_RECORDS_DUMPED = '\n'.join([
                         "expect": 0.11,
                         "frame": 1,
                         "identicalCount": 6,
+                        "percentIdentical": None,
                         "positiveCount": 7,
                         "query": "YSCFT-NSEK",
                         "query_end": 276,
@@ -220,6 +249,7 @@ DIAMOND_RECORDS_DUMPED = '\n'.join([
                         "expect": 0.0,
                         "frame": 1,
                         "identicalCount": 7,
+                        "percentIdentical": None,
                         "positiveCount": 8,
                         "query": "GEPFSVYG",
                         "query_end": 306,
@@ -240,6 +270,7 @@ DIAMOND_RECORDS_DUMPED = '\n'.join([
                         "expect": 7e-143,
                         "frame": 1,
                         "identicalCount": 8,
+                        "percentIdentical": None,
                         "positiveCount": 9,
                         "query": "PFSVYGRF",
                         "query_end": 306,
@@ -265,6 +296,7 @@ DIAMOND_RECORDS_DUMPED = '\n'.join([
                         "expect": 0.008,
                         "frame": 1,
                         "identicalCount": 9,
+                        "percentIdentical": None,
                         "positiveCount": 10,
                         "query": "PKELHGLI",
                         "query_end": 118,
@@ -285,6 +317,7 @@ DIAMOND_RECORDS_DUMPED = '\n'.join([
                         "expect": 0.009,
                         "frame": 1,
                         "identicalCount": 10,
+                        "percentIdentical": None,
                         "positiveCount": 11,
                         "query": "CRPTF",
                         "query_end": 293,
@@ -303,6 +336,252 @@ DIAMOND_RECORDS_DUMPED = '\n'.join([
 ]) + '\n'
 
 
+DIAMOND_RECORDS_WITH_PERCENT_IDENTICAL_DUMPED = '\n'.join([
+    dumps({
+        "application": "DIAMOND",
+        "reference": ("Buchfink et al., Fast and Sensitive "
+                      "Protein Alignment using DIAMOND, Nature Methods, "
+                      "12, 59-60 (2015)"),
+        "task": "blastx",
+        "version": "v0.8.23"
+    }, sort_keys=True),
+    dumps({
+        "alignments": [
+            {
+                "hsps": [
+                    {
+                        "bits": 29.6,
+                        "btop": "4",
+                        "expect": 0.003,
+                        "frame": 1,
+                        "identicalCount": None,
+                        "percentIdentical": 5.0,
+                        "positiveCount": None,
+                        "query": "EFII",
+                        "query_end": 295,
+                        "query_start": 178,
+                        "sbjct": "SSSEV",
+                        "sbjct_end": 285,
+                        "sbjct_start": 175
+                    }
+                ],
+                "length": 295,
+                "title": "INSV"
+            },
+            {
+                "hsps": [
+                    {
+                        "bits": 28.1,
+                        "btop": "3",
+                        "expect": 0.008,
+                        "frame": 1,
+                        "identicalCount": None,
+                        "percentIdentical": 10.0,
+                        "positiveCount": None,
+                        "query": "KLL",
+                        "query_end": 37,
+                        "query_start": 7,
+                        "sbjct": "ITRV",
+                        "sbjct_end": 39,
+                        "sbjct_start": 9
+                    }
+                ],
+                "length": 300,
+                "title": "CASV"
+            },
+            {
+                "hsps": [
+                    {
+                        "bits": 28.1,
+                        "btop": "5",
+                        "expect": 0.009,
+                        "frame": 1,
+                        "identicalCount": None,
+                        "percentIdentical": 15.0,
+                        "positiveCount": None,
+                        "query": "IKSKL",
+                        "query_end": 35,
+                        "query_start": 7,
+                        "sbjct": "EETSR",
+                        "sbjct_end": 37,
+                        "sbjct_start": 9
+                    },
+                    {
+                        "bits": 23.5,
+                        "btop": "6",
+                        "expect": 0.21,
+                        "frame": 1,
+                        "identicalCount": None,
+                        "percentIdentical": 20.0,
+                        "positiveCount": None,
+                        "query": "TIMSVV",
+                        "query_end": 240,
+                        "query_start": 177,
+                        "sbjct": "DDMV",
+                        "sbjct_end": 235,
+                        "sbjct_start": 179
+                    }
+                ],
+                "length": 293,
+                "title": "GoldenGate"
+            },
+            {
+                "hsps": [
+                    {
+                        "bits": 25.0,
+                        "btop": "6",
+                        "expect": 0.084,
+                        "frame": 1,
+                        "identicalCount": None,
+                        "percentIdentical": 25.0,
+                        "positiveCount": None,
+                        "query": "LHVNYL",
+                        "query_end": 203,
+                        "query_start": 1,
+                        "sbjct": "DEELKA",
+                        "sbjct_end": 210,
+                        "sbjct_start": 2
+                    },
+                    {
+                        "bits": 18.5,
+                        "btop": "9",
+                        "expect": 9.1,
+                        "frame": 1,
+                        "identicalCount": None,
+                        "percentIdentical": 30.0,
+                        "positiveCount": None,
+                        "query": "SEIICEVLK",
+                        "query_end": 257,
+                        "query_start": 226,
+                        "sbjct": "VETVAQ",
+                        "sbjct_end": 45,
+                        "sbjct_start": 20
+                    }
+                ],
+                "length": 290,
+                "title": "InfluenzaC"
+            },
+            {
+                "hsps": [
+                    {
+                        "bits": 24.6,
+                        "btop": "10",
+                        "expect": 0.11,
+                        "frame": 1,
+                        "identicalCount": None,
+                        "percentIdentical": 35.0,
+                        "positiveCount": None,
+                        "query": "YSCFT-NSEK",
+                        "query_end": 276,
+                        "query_start": 176,
+                        "sbjct": "LGKRMFC",
+                        "sbjct_end": 243,
+                        "sbjct_start": 152
+                    }
+                ],
+                "length": 270,
+                "title": "FERV"
+            }
+        ],
+        "query": "ACC94"
+    }, sort_keys=True),
+    dumps({
+        "alignments": [
+            {
+                "hsps": [
+                    {
+                        "bits": 634.0,
+                        "btop": "8",
+                        "expect": 0.0,
+                        "frame": 1,
+                        "identicalCount": None,
+                        "percentIdentical": 40.0,
+                        "positiveCount": None,
+                        "query": "GEPFSVYG",
+                        "query_end": 306,
+                        "query_start": 1,
+                        "sbjct": "NIYGEP",
+                        "sbjct_end": 306,
+                        "sbjct_start": 1
+                    }
+                ],
+                "length": 306,
+                "title": "AKAV"
+            },
+            {
+                "hsps": [
+                    {
+                        "bits": 401.0,
+                        "btop": "8",
+                        "expect": 7e-143,
+                        "frame": 1,
+                        "identicalCount": None,
+                        "percentIdentical": 45.0,
+                        "positiveCount": None,
+                        "query": "PFSVYGRF",
+                        "query_end": 306,
+                        "query_start": 1,
+                        "sbjct": "GEPMS",
+                        "sbjct_end": 294,
+                        "sbjct_start": 1
+                    }
+                ],
+                "length": 294,
+                "title": "WYOV"
+            }
+        ],
+        "query": "AKAV"
+    }, sort_keys=True),
+    dumps({
+        "alignments": [
+            {
+                "hsps": [
+                    {
+                        "bits": 28.1,
+                        "btop": "8",
+                        "expect": 0.008,
+                        "frame": 1,
+                        "identicalCount": None,
+                        "percentIdentical": 50.0,
+                        "positiveCount": None,
+                        "query": "PKELHGLI",
+                        "query_end": 118,
+                        "query_start": 14,
+                        "sbjct": "SLKSKE",
+                        "sbjct_end": 131,
+                        "sbjct_start": 15
+                    }
+                ],
+                "length": 307,
+                "title": "TAIV"
+            },
+            {
+                "hsps": [
+                    {
+                        "bits": 28.1,
+                        "btop": "5",
+                        "expect": 0.009,
+                        "frame": 1,
+                        "identicalCount": None,
+                        "percentIdentical": 55.0,
+                        "positiveCount": None,
+                        "query": "CRPTF",
+                        "query_end": 293,
+                        "query_start": 4,
+                        "sbjct": "EFVFIY",
+                        "sbjct_end": 342,
+                        "sbjct_start": 6
+                    }
+                ],
+                "length": 343,
+                "title": "SouthBay"
+            }
+        ],
+        "query": "BHAV"
+    }, sort_keys=True)
+]) + '\n'
+
+
 class TestDiamondTabularFormatReader(TestCase):
     """
     Test the DiamondTabularFormatReader class.
@@ -345,7 +624,7 @@ class TestDiamondTabularFormatReader(TestCase):
     def testDiamondInputWithoutNidentOrPositives(self):
         """
         Test conversion of a chunk of DIAMOND output that does not contain the
-        nident or positives fields.
+        nident, positives, or pident fields.
         """
         mockOpener = mockOpen(
             read_data=DIAMOND_RECORDS_WITHOUT_NIDENT_AND_POSITIVE)
@@ -356,6 +635,8 @@ class TestDiamondTabularFormatReader(TestCase):
                 for alignment in record['alignments']:
                     for hsp in alignment['hsps']:
                         self.assertIs(None, hsp['identicalCount'])
+                        self.assertIs(None, hsp['positiveCount'])
+                        self.assertIs(None, hsp['percentIdentical'])
 
     def testSaveAsJSON(self):
         """
@@ -369,6 +650,20 @@ class TestDiamondTabularFormatReader(TestCase):
             self.maxDiff = None
             self.assertEqual(DIAMOND_RECORDS_DUMPED, fp.getvalue())
 
+    def testSaveAsJSONWithPercentIdentical(self):
+        """
+        A DiamondTabularFormatReader must be able to save itself as JSON
+        when the percentIdentical field is present.
+        """
+        mockOpener = mockOpen(read_data=DIAMOND_RECORDS_WITH_PERCENT_IDENTICAL)
+        with patch.object(builtins, 'open', mockOpener):
+            reader = DiamondTabularFormatReader('file.txt')
+            fp = StringIO()
+            reader.saveAsJSON(fp)
+            self.maxDiff = None
+            self.assertEqual(DIAMOND_RECORDS_WITH_PERCENT_IDENTICAL_DUMPED,
+                             fp.getvalue())
+
     def testSaveAsJSONBzip2(self):
         """
         A DiamondTabularFormatReader must be able to save itself as bzip2'd
@@ -385,6 +680,23 @@ class TestDiamondTabularFormatReader(TestCase):
                 compress(DIAMOND_RECORDS_DUMPED.encode('UTF-8')),
                 data.getvalue())
 
+    def testSaveAsJSONBzip2WithPercentIdentical(self):
+        """
+        A DiamondTabularFormatReader must be able to save itself as bzip2'd
+        JSON when the percentIdentical field is present.
+        """
+        mockOpener = mockOpen(read_data=DIAMOND_RECORDS_WITH_PERCENT_IDENTICAL)
+        with patch.object(builtins, 'open', mockOpener):
+            reader = DiamondTabularFormatReader('file.txt')
+            data = BytesIO()
+            fp = bz2file.BZ2File(data, 'w')
+            reader.saveAsJSON(fp, writeBytes=True)
+            fp.close()
+            self.assertEqual(
+                compress(DIAMOND_RECORDS_WITH_PERCENT_IDENTICAL_DUMPED.encode(
+                    'UTF-8')),
+                data.getvalue())
+
     def testSpacesMustBePreserved(self):
         """
         If there are spaces in the query title or subject titles, the spaces
diff --git a/test/test_hsp.py b/test/test_hsp.py
index 2494ec3..b9b645e 100644
--- a/test/test_hsp.py
+++ b/test/test_hsp.py
@@ -64,7 +64,7 @@ class TestHSP(TestCase):
                   subjectStart=5, subjectEnd=6,
                   readMatchedSequence='aaa', subjectMatchedSequence='ccc',
                   readFrame=7, subjectFrame=8, identicalCount=9,
-                  positiveCount=10)
+                  positiveCount=10, percentIdentical=99.3)
 
         self.assertEqual(
             {
@@ -81,6 +81,7 @@ class TestHSP(TestCase):
                 'positiveCount': 10,
                 'readMatchedSequence': 'aaa',
                 'subjectMatchedSequence': 'ccc',
+                'percentIdentical': 99.3,
             },
             hsp.toDict())
 
@@ -146,7 +147,7 @@ class TestLSP(TestCase):
                   subjectStart=5, subjectEnd=6,
                   readMatchedSequence='aaa', subjectMatchedSequence='ccc',
                   readFrame=7, subjectFrame=8, identicalCount=9,
-                  positiveCount=10)
+                  positiveCount=10, percentIdentical=99.3)
 
         self.assertEqual(
             {
@@ -163,5 +164,6 @@ class TestLSP(TestCase):
                 'positiveCount': 10,
                 'readMatchedSequence': 'aaa',
                 'subjectMatchedSequence': 'ccc',
+                'percentIdentical': 99.3,
             },
             lsp.toDict())
diff --git a/test/test_titles.py b/test/test_titles.py
index d170077..53c0a2c 100644
--- a/test/test_titles.py
+++ b/test/test_titles.py
@@ -58,13 +58,13 @@ class TestTitleAlignment(TestCase):
                    subjectStart=5, subjectEnd=6,
                    readMatchedSequence='aaa', subjectMatchedSequence='ccc',
                    readFrame=7, subjectFrame=8, identicalCount=9,
-                   positiveCount=10)
+                   positiveCount=10, percentIdentical=31.3)
         hsp2 = HSP(10, readStart=11, readEnd=12,
                    readStartInSubject=13, readEndInSubject=14,
                    subjectStart=15, subjectEnd=16,
                    readMatchedSequence='ggg', subjectMatchedSequence='ttt',
                    readFrame=17, subjectFrame=18, identicalCount=19,
-                   positiveCount=20)
+                   positiveCount=20, percentIdentical=32.3)
         titleAlignment = TitleAlignment(read, [hsp1, hsp2])
 
         self.assertEqual(
@@ -84,6 +84,7 @@ class TestTitleAlignment(TestCase):
                         'positiveCount': 10,
                         'readMatchedSequence': 'aaa',
                         'subjectMatchedSequence': 'ccc',
+                        'percentIdentical': 31.3,
                     },
                     {
                         'score': 10,
@@ -99,6 +100,7 @@ class TestTitleAlignment(TestCase):
                         'positiveCount': 20,
                         'readMatchedSequence': 'ggg',
                         'subjectMatchedSequence': 'ttt',
+                        'percentIdentical': 32.3,
                     },
                 ],
                 'read': {
@@ -899,13 +901,13 @@ class TestTitleAlignments(WarningTestMixin, TestCase):
                    subjectStart=5, subjectEnd=6,
                    readMatchedSequence='aaa', subjectMatchedSequence='ccc',
                    readFrame=7, subjectFrame=8, identicalCount=9,
-                   positiveCount=10)
+                   positiveCount=10, percentIdentical=17.9)
         hsp2 = HSP(10, readStart=11, readEnd=12,
                    readStartInSubject=13, readEndInSubject=14,
                    subjectStart=15, subjectEnd=16,
                    readMatchedSequence='ggg', subjectMatchedSequence='ttt',
                    readFrame=17, subjectFrame=18, identicalCount=19,
-                   positiveCount=20)
+                   positiveCount=20, percentIdentical=27.9)
         titleAlignment = TitleAlignment(read, [hsp1, hsp2])
         titleAlignments = TitleAlignments('subject title', 10)
         titleAlignments.addAlignment(titleAlignment)
@@ -931,6 +933,7 @@ class TestTitleAlignments(WarningTestMixin, TestCase):
                                 'positiveCount': 10,
                                 'readMatchedSequence': 'aaa',
                                 'subjectMatchedSequence': 'ccc',
+                                'percentIdentical': 17.9,
                             },
                             {
                                 'score': 10,
@@ -946,6 +949,7 @@ class TestTitleAlignments(WarningTestMixin, TestCase):
                                 'positiveCount': 20,
                                 'readMatchedSequence': 'ggg',
                                 'subjectMatchedSequence': 'ttt',
+                                'percentIdentical': 27.9,
                             },
                         ],
                         'read': {