Finding homologous regions in multiple whole genomes

Question

Right now I have six genomes that I want to compare and identify homologous regions in the genomes. I have run nucmer, show-coords and obtained the output files. An example is shown below with Genome 1 vs Genome 2. More files go through Genome 1 vs Genome 3 until all are compared.

I was trying to use bedtools merge after concatenating the multiple output files from show-coords command to find the homologous regions. However, I am not quite sure on the command to use. Any help is appreciated.

Edit: added commands. I just want to know what regions are homologous across the six genomes (99/100% IDY). I don't know what I am supposed to merge to find this out. I was reading and I believe bedtools should be but I don't know. I am usually a wetlab but I try to do a bit of bioinformatics.

 # Loop through all pairs of genomes
for i in {1..6}; do
    for j in $(seq $((i + 1)) 6); do
        genome1="genome${i}.fasta"
        genome2="genome${j}.fasta"

        output_prefix="output_${i}_${j}"

        # Run NUCmer
        nucmer -p "${output_prefix}" "$genome1" "$genome2"

        # Show coordinates
        show-coords -r -l -c -T "${output_prefix}.delta" > "${output_prefix}_coords.txt"
    done
done

After using the above I then concatenate the files with

cat output_*_coords.txt > combined_coords.txt

This gives me something that looks like this below.

[S1]    [E1]    [S2]    [E2]    [LEN 1] [LEN 2] [% IDY] [LEN R] [LEN Q] [COV R] [COV Q] [TAGS]
1       14347   1       14347   14347   14347   99.02   5064019 5004926 0.28    0.29    Genome1      Genome2
14474   32175   14350   32051   17702   17702   99.33   5064019 5004926 0.35    0.35    Genome1      Genome2
16541   22042   4698474 4692975 5502    5500    99.64   5064019 5004926 0.11    0.11    Genome1      Genome2
16548   22042   3728461 3722968 5495    5494    99.84   5064019 5004926 0.11    0.11    Genome1      Genome2
16583   22041   3904238 3898781 5459    5458    99.65   5064019 5004926 0.11    0.11    Genome1      Genome2
18727   22097   4426876 4423506 3371    3371    99.73   5064019 5004926 0.07    0.07    Genome1      Genome2

Answer 1

According to show cords the format you are getting contains columns which you can use to obtain BED format with relevant information.
You need following columns:

TAGS_genome1 
[S1]    Start of the alignment region in the reference sequence.
[E1]    End of the alignment region in the reference sequence.
TAGS_genome2
[% IDY]

I would start with not merging all coordinates but getting bed files for each pairwise genome_a__vs__genome_b.bed etc. and start from there, taking a subset of 3 genomes first (say a, b,c). You can find the regions shared by all three bedtools intersect:

bedtools intersect -a a_vs_b.bed -b a_vs_c.bed

The issue is that you will have many to many relationships, since region_1 from b can have several similar regions in b, and each one of these (== regions in b) also can be similar to a number of regions in a. Depending on what you are looking for and how dis-similar the genomes are you may end up with reporting repetitive elements or if you filter hard missing i.e. exons.