Various updates #421

nh13 · 2018-07-28T02:40:23Z

@tfenne here are some various updates I have been accumulating.

codecov-io · 2018-07-28T02:51:30Z

Codecov Report

Merging #421 into master will decrease coverage by 0.01%.
The diff coverage is 96.15%.

@@            Coverage Diff             @@
##           master     #421      +/-   ##
==========================================
- Coverage   96.27%   96.26%   -0.02%     
==========================================
  Files          90       90              
  Lines        5099     5109      +10     
  Branches      659      652       -7     
==========================================
+ Hits         4909     4918       +9     
- Misses        190      191       +1

Impacted Files	Coverage Δ
...ava/com/fulcrumgenomics/util/PickLongIndices.scala	`94.49% <100%> (ø)`	⬆️
...ain/scala/com/fulcrumgenomics/util/Sequences.scala	`100% <100%> (ø)`	⬆️
...c/main/scala/com/fulcrumgenomics/util/Metric.scala	`93.75% <100%> (+0.2%)`	⬆️
...cala/com/fulcrumgenomics/util/ProgressLogger.scala	`50% <66.66%> (-16.67%)`	⬇️

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 431b695...3397e55. Read the comment docs.

tfenne

Sorry for being so picky :/

tfenne · 2018-07-30T13:40:06Z

src/main/scala/com/fulcrumgenomics/util/ProgressLogger.scala

+  def record(rec: SamRecord): Boolean = record(rec.asSam)
+
+  /** Logs the last record if it wasn't already logged. */
+  def logLast(): Boolean = super.log()


It's really ugly that we have both log() and log(String... message) in the superclass. So when you refer to log() it's not very clear what's being called :/ Can you maybe add a trailing comment here something like:

def logLast(): Boolean = super.log() // Calls the super's log() method, not log(message: String*)

tfenne · 2018-07-30T13:46:20Z

src/main/scala/com/fulcrumgenomics/util/Sequences.scala

+  }
+
+
+  case class OffsetAndLength(offset: Int, length: Int)


Should have at least a single line scaladoc if it's a public class.

tfenne · 2018-07-30T13:50:05Z

src/main/scala/com/fulcrumgenomics/util/Sequences.scala

+    * @return the offset and length of the longest homopolymer
+    */
+  def maxPolyX(s: String) : OffsetAndLength = {
+    val bs = s.getBytes


Doing this vs. using charAt() causes the creation of a whole new byte array every time this function is called.

tfenne · 2018-07-30T13:51:58Z

src/main/scala/com/fulcrumgenomics/util/Sequences.scala

+    val bs = s.getBytes
+    var (bestStart, bestLength) = (0,0)
+    var start = 0
+    while (start < bs.length) {


Why switch to using a while loop instead of the forloop that was used previously? It automatically makes me think you're doing something other than incrementing by 1 each iteration.

tfenne · 2018-07-30T13:52:29Z

src/main/scala/com/fulcrumgenomics/util/Sequences.scala

+    var (bestStart, bestLength) = (0,0)
+    var start = 0
+    while (start < bs.length) {
+      val firstBase = s.charAt(start).toByte


Odd to use s.charAt(start) here and bs(i) below.

tfenne · 2018-07-30T14:00:04Z

src/main/scala/com/fulcrumgenomics/util/Sequences.scala

+    forloop(0, bs.length-1) { start =>
+      val (b1, b2) = (bs(start), bs(start+1))
+      var i = start
+      while (i < s.length-1 && SequenceUtil.basesEqual(b1, bs(i)) && SequenceUtil.basesEqual(b2, bs(i+1))) i += 2


Two thoughts:

If you upper case like you're doing on 105 you don't really need to use SequenceUtil.basesEqual() here.

If you continue to use SequenceUtil.basesEqual() I'd be tempted to define and use the following within this class:

/** Checks to see if the two characters are identical, ignoring case. */ @inline final private def same(a: Char, b: Char): Boolean = { Character.toUpperCase(a) == Character.toUpperCase(b) }

This has the advantage that a) it'll work for all characters so if you try to use it on proteins it won't fail, and also will just make the usage more compact and easier to see what's going on.

Done, but "I did it my way"

tfenne · 2018-07-30T14:03:45Z

src/main/scala/com/fulcrumgenomics/util/Sequences.scala

+  /** Returns the sequence that is the complement of the provided sequence. */
+  def complement(s: String): String = {
+    val bs = s.getBytes
+    for (i <- 0 until bs.length) bs(i) = SequenceUtil.complement(bs(i))


I think forloop(from=0, until=bs.length) will be much faster here.

It is, fixed.

tfenne · 2018-07-30T14:05:23Z

src/test/scala/com/fulcrumgenomics/util/SequencesTest.scala

+  }
+
+  it should "return zero when passed a zero length string" in {
+    Sequences.maxPolyX("") shouldBe OffsetAndLength(0, 0)


Why would you want this behavior? Also, I don't think it's reasonable to deprecate one method in favor of another if the behavior is different.

Sometimes I have case classes with an empty sequence, so I don't want it to blow up when calling maxPolyX:

things.map(maxPolyX(_.sequence))

The longest homopolymer in an empty sequence is the empty string isn't it?
2. I don't think you should be using this method to catch the case (i.e. throw an exception) when you have an empty string. I think it's wrong, so I'd like to fix it. I really don't want to have to write guards everywhere in my code to ignore feeding empty strings to longestHomopolymer. One could argue that asking for the longest homopolymer implies that there is at least one homopolymer, but I think that's just an argument for argument's sake.

Any objections to fixing this? Do you have code that will break (i.e. relies on the exception being thrown here)?

tfenne · 2018-07-30T14:09:34Z

src/main/scala/com/fulcrumgenomics/util/Sequences.scala

+    * @param s a DNA or RNA sequence
+    * @return the offset and length of the longest homopolymer
+    */
+  def maxPolyX(s: String) : OffsetAndLength = {


Two thoughts on naming (apply to maxDinux as well):

I think longest is clearer than max. Longest can only really be interpreted one way. Max could be interpreted to mean "longest single", or "most abundant" such that in a sequence like ACACGTGTGTACAC you might say AC is the "max dinuc" because there is more of it.

Homopolymer is a well defined term whereas polyX isn't really.

Agreed, fixed.

tfenne · 2018-07-30T14:11:34Z

src/test/scala/com/fulcrumgenomics/util/SequencesTest.scala

+    Sequences.complement("AAA") shouldBe "TTT"
+    Sequences.complement("ACAC") shouldBe "TGTG"
+    Sequences.complement("") shouldBe ""
+    Sequences.complement("AACCGGTGTG") shouldBe "TTGGCCACAC"
  }


I know it's only a call to a htsjdk function, but could you add a few basic tests for reverseComplement() please? That way if the library function gets broken/changed unexpectedly we'll catch it, or if we decide to reimplement we won't be caught short.

nh13

@tfenne I made all your suggested fixes with the exception of the longestHomopolymer method failing when an empty string is given. Lets slack chat if you object to my reasoning to allow an empty string.

nh13 · 2018-08-11T22:44:55Z

src/test/scala/com/fulcrumgenomics/util/SequencesTest.scala

+  }
+
+  it should "return zero when passed a zero length string" in {
+    Sequences.maxPolyX("") shouldBe OffsetAndLength(0, 0)


Sometimes I have case classes with an empty sequence, so I don't want it to blow up when calling maxPolyX:

things.map(maxPolyX(_.sequence))

The longest homopolymer in an empty sequence is the empty string isn't it?
2. I don't think you should be using this method to catch the case (i.e. throw an exception) when you have an empty string. I think it's wrong, so I'd like to fix it. I really don't want to have to write guards everywhere in my code to ignore feeding empty strings to longestHomopolymer. One could argue that asking for the longest homopolymer implies that there is at least one homopolymer, but I think that's just an argument for argument's sake.

Any objections to fixing this? Do you have code that will break (i.e. relies on the exception being thrown here)?

nh13 · 2018-08-11T22:45:02Z

src/test/scala/com/fulcrumgenomics/util/SequencesTest.scala

+    Sequences.complement("AAA") shouldBe "TTT"
+    Sequences.complement("ACAC") shouldBe "TGTG"
+    Sequences.complement("") shouldBe ""
+    Sequences.complement("AACCGGTGTG") shouldBe "TTGGCCACAC"
  }


nh13 · 2018-08-11T22:45:09Z

src/main/scala/com/fulcrumgenomics/util/ProgressLogger.scala

+  def record(rec: SamRecord): Boolean = record(rec.asSam)
+
+  /** Logs the last record if it wasn't already logged. */
+  def logLast(): Boolean = super.log()


nh13 · 2018-08-11T22:45:13Z

src/main/scala/com/fulcrumgenomics/util/Sequences.scala

+  }
+
+
+  case class OffsetAndLength(offset: Int, length: Int)


nh13 · 2018-08-11T22:45:18Z

src/main/scala/com/fulcrumgenomics/util/Sequences.scala

+    * @param s a DNA or RNA sequence
+    * @return the offset and length of the longest homopolymer
+    */
+  def maxPolyX(s: String) : OffsetAndLength = {


Agreed, fixed.

nh13 · 2018-08-11T22:45:30Z

src/main/scala/com/fulcrumgenomics/util/Sequences.scala

+    var (bestStart, bestLength) = (0,0)
+    var start = 0
+    while (start < bs.length) {
+      val firstBase = s.charAt(start).toByte


nh13 · 2018-08-11T22:45:33Z

src/main/scala/com/fulcrumgenomics/util/Sequences.scala

+    * @return the offset and length of the longest dinucleotide sequence
+    */
+  def maxDinuc(s: String) : OffsetAndLength = {
+    val bs = s.toUpperCase.getBytes


nh13 · 2018-08-11T22:45:45Z

src/main/scala/com/fulcrumgenomics/util/Sequences.scala

+    forloop(0, bs.length-1) { start =>
+      val (b1, b2) = (bs(start), bs(start+1))
+      var i = start
+      while (i < s.length-1 && SequenceUtil.basesEqual(b1, bs(i)) && SequenceUtil.basesEqual(b2, bs(i+1))) i += 2


Done, but "I did it my way"

nh13 · 2018-08-11T22:45:52Z

src/main/scala/com/fulcrumgenomics/util/Sequences.scala

+  /** Returns the sequence that is the complement of the provided sequence. */
+  def complement(s: String): String = {
+    val bs = s.getBytes
+    for (i <- 0 until bs.length) bs(i) = SequenceUtil.complement(bs(i))


It is, fixed.

nh13 · 2018-08-11T22:46:09Z

src/main/scala/com/fulcrumgenomics/util/Sequences.scala

+  }
+
+  /** Reverse complements a string of bases. */
+  def reverseComplement(s: String): String = SequenceUtil.reverseComplement(s)


Feeling good, you? Fixed.

This is useful if you want to log the final record (or count).

In some cases, the metric file would not be fully written.

Approval over Slack

nh13 assigned tfenne Jul 28, 2018

nh13 requested a review from tfenne July 28, 2018 02:40

nh13 added the enhancement label Jul 29, 2018

tfenne previously requested changes Jul 30, 2018

View reviewed changes

tfenne assigned nh13 and unassigned tfenne Jul 30, 2018

nh13 commented Aug 11, 2018

View reviewed changes

nh13 assigned tfenne and unassigned nh13 Aug 11, 2018

nh13 added 3 commits August 13, 2018 14:57

ProgressLogger.record should return if logging was triggered

963f302

Add a method to log the last record in ProgressLogger.

5aa574b

This is useful if you want to log the final record (or count).

Bug fix: Metric.write was not closing its writer.

1cd2ddd

In some cases, the metric file would not be fully written.

nh13 force-pushed the nh_various_updates branch 2 times, most recently from 27ae815 to 81dea48 Compare August 13, 2018 21:59

Adding a few useful methods to Sequences.

3397e55

nh13 force-pushed the nh_various_updates branch from 81dea48 to 3397e55 Compare August 13, 2018 22:06

nh13 merged commit f2d5d5e into master Aug 13, 2018

nh13 deleted the nh_various_updates branch August 13, 2018 22:36

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Various updates #421

Various updates #421

nh13 commented Jul 28, 2018

codecov-io commented Jul 28, 2018 •

edited

Loading

tfenne left a comment

tfenne Jul 30, 2018

nh13 Aug 11, 2018

tfenne Jul 30, 2018

nh13 Aug 11, 2018

tfenne Jul 30, 2018

nh13 Aug 11, 2018

tfenne Jul 30, 2018

nh13 Aug 11, 2018

tfenne Jul 30, 2018

nh13 Aug 11, 2018

tfenne Jul 30, 2018

nh13 Aug 11, 2018

tfenne Jul 30, 2018

nh13 Aug 11, 2018

tfenne Jul 30, 2018

nh13 Aug 11, 2018

tfenne Jul 30, 2018

nh13 Aug 11, 2018

tfenne Jul 30, 2018

nh13 Aug 11, 2018

nh13 left a comment

nh13 Aug 11, 2018

nh13 Aug 11, 2018

nh13 Aug 11, 2018

nh13 Aug 11, 2018

nh13 Aug 11, 2018

nh13 Aug 11, 2018

nh13 Aug 11, 2018

nh13 Aug 11, 2018

nh13 Aug 11, 2018

nh13 Aug 11, 2018

Various updates #421

Various updates #421

Conversation

nh13 commented Jul 28, 2018

codecov-io commented Jul 28, 2018 • edited Loading

Codecov Report

tfenne left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

nh13 left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

codecov-io commented Jul 28, 2018 •

edited

Loading