import pybedtools a = pybedtools.BedTool('H3K4me1_enhancers.bed') b = pybedtools.BedTool('H3K27ac_enhancers.bed') c = pybedtools.BedTool('GRO-seq_filtered_peaks.bed') d = (a+b+c) e = (a+b-c) f = (c-a-b) d.saveas('All_enhancers.bed') e.saveas('Histone_only_enhancers.bed') f.saveas('GRO_only_enhancers.bed')