code: def coverage(self):
          
          variant_start_index, variant_end_index = self.variant_indices()
          n_nucleotides = len(self)
          coverage_array = np.zeros(n_nucleotides, dtype="int32")
          for read in self.reads:
              coverage_array[
                  max(0, variant_start_index - len(read.prefix)):
                  min(n_nucleotides, variant_end_index + len(read.suffix))] += 1
          return coverage_array
docs: Returns NumPy array indicating number of reads covering each
          nucleotides of this sequence.