If it's okay to use BEDOPS, you could use a one-liner of bedops --partition, bedmap and awk:

$ bedops --partition regions.bed |\
  bedmap --echo --echo-map-id-uniq --delim '\t' - regions.bed |\
  awk '{
          n = split($4, a, ";");
          max = a[1];
          for(i = 2; i <= n; i++)
          {
              if (a[i] > max)
              {
                  max = a[i];
              }
          }
          print $1"\t"$2"\t"$3"\t"max;
      }' - > answer.bed

Make sure the input regions are tab-delimited and sorted with sort-bed.

This works by partitioning the input into disjoint regions however many overlaps there are (bedops), mapping the disjoint regions to the original input to get a string of overlapping ids (bedmap), and using awk to find the maximum value in that string. Then awk prints the disjoint region and its max id.

By using standard input/output streams and using a linear sweep for the line maximum, this should run fast and have low memory overhead.

Here's what things look like at each step of the pipeline.

First, we partition the input regions into disjoint elements:

$ bedops --partition regions.bed
chr1    0      50
chr1    50     100
chr1    100    150
chr1    150    200

You can throw as many overlaps into --partition as you like.

Then we map the partitioned set back to the original input regions, and get the list of overlapping ids:

$ bedops --partition regions.bed |\
  bedmap --echo --echo-map-id-uniq --delim '\t' - regions.bed
chr1    0      50     2
chr1    50     100    1;2
chr1    100    150    1;3
chr1    150    200    3

Then awk parses each of the id-strings in the fourth column to get the maximum value, printing the result:

$ bedops --partition regions.bed |\
  bedmap --echo --echo-map-id-uniq --delim '\t' - regions.bed |\
  awk '{ 
          n = split($4, a, ";");
          max = a[1];
          for(i = 2; i <= n; i++)
          {
              if (a[i] > max)
              {
                  max = a[i];
              }
          }
          print $1"\t"$2"\t"$3"\t"max;
       }' -
chr1    0      50     2
chr1    50     100    2
chr1    100    150    3
chr1    150    200    3

I believe you cannot solve this without writing your own algorithm. Here's my first attempt:

cat a.bed

chr1    0       100     2
chr1    50      150     1
chr1    100     200     3

(cut -f2 a.bed; cut -f3 a.bed) |\
    sort -nu |\
    perl -ne '
        BEGIN{$i=0} chomp;
        push @F, $_;
        if($i++){print "chr1\t$F[$i-2]\t$F[$i-1]\n"}
        ' > b.bed

bedtools intersect -a a.bed -b b.bed | sort -k2n -k3n -k4nr | perl -lane 'print unless $h{$F[0,1,2]}++'

Output:

chr1    0       50      2
chr1    50      100     2
chr1    100     150     3
chr1    150     200     3

You can try doing this:

1. Use clusterBed to cluster overlapping intervals

2. For each cluster get the elementary intervals. Basically just put all start,end coordinates into a list, make the list unique, then sort the list. So something like this:

   for cluster in clusters:
    coords = set()
    for line in cluster:
        cols = line.strip().split()
        ref = cols[0]
        start, end, cov = map(int, cols[1:])
        coords.add(start)
        coords.add(end)
    coords = list(coords)
    coords.sort()
    elementary_intervals = []
    for i in range(len(coords) - 1):
        elementary_intervals.append((coords[i],coords[i+1]))
   

3. For each elementary interval, find overlaps with cluster members and find the maximum coverage (or whatever other operations). So something like this:

   for e_interval in e_intervals:
    covs = []
    for member in cluster:
        if overlap(e_interval,member):
            covs.append(memberCov)
    max(covs) or sum(covs) or whatever