同样的输入文件多次运行得到的结果不一致

[Threeman8]

您好！非常感谢您的工作！

我用同样的输入文件，hic、hifi和genome，其他参数默认，运行了两次，得出来的结果不一致。我将两次得到的数量最多的CL进行了blastn，它们相似性也很低。
第一次最多的CL是：
`

CL8 N=3;L=273;center=TR250;weight=0.36;ratio=0.002%
AAGATGGTGCCACCTCTAAACGACAACCTAAACGGACCAAAACGTGACGAAAATGGTGCC
ACTCTCTTCACAACGGTTTGTGGAATATGCGACTCTTTCATCGCACAAAGAGGAGGAAAC
GGGACGAAGATGATGCCACATCTAAACGACACCCTAAACGGACCAAAACGTGACAAAAAT
GGTACCACTCTCTTCACAACGGTTTGTGGAATATTCGACTATTTCATCGCACAAAGAGGA
GGAAACGAGAGCCCAAAGTGGAGGAAACGGGCG
CL6 N=4;L=147;center=TR332;weight=0.26;ratio=0.002%
ATGGTGCTACCTTTAAACGACAACCTAAACGGACCAAAATGTGACGAAAATGGTGCACTC
TCTTCACAACGGTTTGTGGACTATGCGACTCTTTCATCGCACAAAGAGGAGGAAACGGGA
GCCCAAAGTGGAGGAAACGGGACGAAG
`

第二次最多的CL是：
`

CL6 N=4;L=147;center=TR227;weight=0.26;ratio=0.002%
ATGGTGCTACCTTTAAACGACAACCTAAACGGACCAAAATGTGACGAAAATGGTGCACTC
TCTTCACAACGGTTTGTGGACTATGCGACTCTTTCATCGCACAAAGAGGAGGAAACGGGA
GCCCAAAGTGGAGGAAACGGGACGAAG

CL10 N=2;L=297;center=TR545;weight=0.21;ratio=0.001%
TGAAGAGAGTGCACCATTTTCATCATGTTTTGGTCCGTTTAGAGTGCCGTTTAAAGGTAG
CACCATATTCGTCCCGTGTCCTCCACTTTGGGCTCCCGTTTGCTCCTCTTTGTGCACGAT
GAAAGAGTCGCATAATCCACAAACCGTTGTGAATAGAGTGCACCATTTTCGTCACGTTTT
GGTCAGTTTAATGTGTCGTTTAAATGTAGCACCATCTTCGTCCTGTTTCCTCCACTTTGG
GGTCCCGTTTGCTCCTCTTTGTGCGGATGAAAGAATTGCATAGTCCATAAACCATTG

CL8 N=2;L=1756;center=TR178;weight=0.74;ratio=0.001%
ACGGGAGCCCAAAGTGGAGGAAATGGGACGAAGATGGTGCTACCTTTAAACGACAACCTA
AACAAACAAAAACGTGACAAAAATGGTGCACTCTCTTCACATCGGTTTGCGGACTATGCA
ACTCTTTTAGCACACACAAAGAGGAACAAACCAGAGCCCAAAGTGGAGGAAACAGGACGA
AGATGGTGCTACCTTTAAACGACACACTAAACGGATCAAAACGTGATGAAAATGCTGCAA
TCTCTTACAATGGTTTGTGGACTATGCTACTCTTTAATCGTGCACAAAGAGGAGCAAACG
GGAGCCCAAAGTGGAGGAAATGGGACGAAGGTGGTGCTACTTTTAATCGAAACACTAAAT
GGACCAAAACGTGATGAAAATGGTGCAACTCTAATCACAATGGTTTTTGGACTATGTGAC
TCTTTCATCGTACATAGGGGAGGAAACGGGATCCCAAAGTGGAGGAAACGGGACAAATAT
GGTGCTACCTTAAAAGACACACTAAACACACCAAAACCTGACGATAATGGTGCATTCTCT
TCACAACGGTTGGTGGACTATGCGACTCTTTCATCGCACACAAAGAGGAGGAAACGGGAC
GAAGATGGTGCTACCTTTAAACGGCAAGCTAAACAGACCAAAATGTGACGAAAATGGTGC
ACTCTATTCACAATGGTTTGTGGAATATGCTACTCTTTCATCGCATAAAGAAAATGAAAC
GGAAGCCCAAAGTGGAGGAAATGGGACGAAGATGGTGCTACATTTAAATGACACACTAAA
CAGACCAAAATATGACAAAAATGGTACCACTCTCTTCACAATGGTTTGAGGAATATGCTA
CTCTTTCATCGCATAAAGAAAATCAAACGGAAGCCCAAAGTGGAGGAAATGGGACGAAGA
TGGTGCTACCTTTAAACGATAAACTAAACAAACCAAAACCTGATGAAAATGGTGCACTCT
CTTCACAACGGTTTGTGGACTATGCAACTCTTTCAGCACGCACAAAGAGGAGTAAACGAG
AGCCCAAATTGGAGGAAACGGGACAAATATAGTGCCACCTTTATACGACCACCTAAATGA
ACCAAAATGTGACAAAAATGGTACCACTCTCTTTACAATGGTTTGTGGAATATGCTACTC
TTTCATCGCATAAAGAGGAGGAAACTGGAGCCCAAAGTGGAGGAAATGGGATGAAGATGG
TGCTACCTTTAAACGACAACCTAAATAAACCAAAACGTGACGAAAATGGTGCACTCTCTT
GACAACGATTGATGCACTATGAAACACTTTTATCACGCACAAAGAATAGCAAACGAGATC
CCAAAATGGAGGAAACACGACGAAGATGGTGCTACCTTTAAACGACACACAAAACGGACC
AAAACGTGATGAAAATGGTGCACACTCTTACAACGGTTTGTGGACTATGCGACTCTTTAA
TCGCGCACAAAGAGGAGGAAACGGGAGCCCAAAGTGGAGGAAACAGGACGATGATGACGC
TACCTTTAAACGACAACCTAAACGAACCAAATGTGACAAAAATGGTGCACTCTCTTCACA
ACGGTTTAAGGACTATGCAACTCTTTCATCGCGCACAAAGAGGAGGAAACAGGAACCCAA
AGTGGAGGAAACGGGACGAAGATGGTGCCACCTTTAAACGACAACCTAAATGAACCAAAA
TGTGACTAAAATGGTACCACTATCTTCACAACGGTTTGTGGAATATGTTACTCTTTCATC
GCATTAAGAAGAGGAA

`

它们在染色体的位置是相似的。

请问这是什么原因呢？
谢谢！

[zhangrengang]

默认是随机取5x的hifi reads（全用可能太慢）去检测串联重复，可能是随机性导致的，可以提高数据量（参数-subsample_x）试试。你这个似乎串联重复的含量很低（容易受随机因素影响），拟南芥能到2%。

[zhangrengang]

你这个N也都很低，说明只有少量几条reads检测到了这个串联重复。

[Threeman8]

好咧，那我指定-subsample_n 试试，稍后回复您。

另，这里的注释是什么意思呀？
CL8 N=2;L=1756;center=TR178;weight=0.74;ratio=0.001%

[zhangrengang]

N是cluster的基序数量，L是基序长度，ratio是累计长度占比，就是在基因组中的大致比例，另两个是聚类相关的。

[Threeman8]

请问ratio是怎么算的呢？
我是这么算：
对cent-output文件夹中的trf.count文件进行：
for i in {4..14} ; do awk -v col=$i '{sum += $col} END {print sum}' trf.count ;echo ${i}; done
拿到count后，乘以motif长度L，然后除以genome.fa的长度。

不知道是哪里出了问题呢？
另，请问script是在/miniconda3/envs/RepCent/bin/centromics这里吗？

[zhangrengang]

ratio是基于reads算的，motif总长度除以reads总长度（motif是从reads中鉴定的）。
什么script？

[Threeman8]

想知道centromics软件运行时，具体是怎么运行得出这些结果的

[zhangrengang]

可以读下代码。。

[Threeman8]

对对，我意思是在哪里可以看到代码

[zhangrengang]

代码就在这个github仓库。

[Threeman8]

感谢感谢！（新手问题有点愚蠢，请见谅）
-subsample_x 80 的结果显示最像着丝粒单元的是CL6和CL3
`>CL3 N=133;L=149;center=TR11362;weight=0.40;ratio=0.005%
GGAGGAAACGGGACGAAGATGGTGCTACCTTTAAACGACACACTAAACGGACCAAAACGT
GACGAAAATGGTGCACTCTCTTCACAACGGTTTGTGGACTATGCGACTCTTTCATCGCGC
ACAAAGAGGAGGAAACGGGAGCCCAAAGT

CL6 N=78;L=148;center=TR2072;weight=0.41;ratio=0.003%
AATATGCGACTCTTTCATCGCACAAAGAGGAGGAAACGGGAGCCCAAAGTGGAGGAAACG
GGACGAAGATGGTGCCACCTCTAAACGACACCCTAAACGGACCAAAACGTGATGAAAATG
GTGCCACTCTCTTCACAACGGTTTGTGG
`
ratio还是挺低的。CL6和CL3的blastn相似度是96%，为啥他们不会聚类成一个呢？

[zhangrengang]

你可以调聚类参数让它们聚成一个。