☰

如何使用 shuf 从随机提取 200 个字符的脚本中进行无替换采样？

tag-icon • tag-icon shell-script shell for random shuf
如何使用 shuf 从随机提取 200 个字符的脚本中进行无替换采样？

我有一个从一组随机字符中提取 200 个字符的脚本：

#!/usr/bin/bash
n=$(stat -c "%s" newfile.txt)
r=$(shuf -i1-"$((n-200+1))" -n1)
< newfile.txt tail -c+"$r" | head -c200

for N in {1..10000}; do bash extraction_200pb.sh; done > output.txt

我知道shuf它非常强大，但我想包含一个不替换的样本。这意味着采样时每抽取200个字符只有一次被选中的机会。

输出应如下所示：

>1     
GAACTCTACCAAAAGGTATGTTGCTTTCACAAAAAGCTGCATTCGATCATGTGTATAATCTAGCAAAACTAGTAGGAGGAGCAAAATACCCCGAAATTGTTGCTGCTCAGGCAATGCACGAATCAAACTACCTAGATCCTAGG
ACTAATAGTGTTTATAATGCCACAAATAGAACTAATGCTTTCGGTCAAACTGGTGAC
>2     
GCCTACCGCATAAAACAGCATCACCGCCACGGCTTCAGGGTATTCTCCAATGGCAAAGGCTCCCATGGTCGCGATGGACATTAAGAGAAATTCAGTAAAGAAATCTCCATTTAGAATACTTTTGAATCCTTCTTTTATCACCG
GAAAACCAACTGGGAGATAGGCCACAATGTACCAACCTACTCGCACCCAATCTGTAA
>3     
GCACGTGTCACCGTCAGCATCGCGGCAGCGGAACGGGTCACCCGGATTGCTGTCGGGACCATCGTTTACGCCGTCATTGTCGTTATCGGGATCGCCCGGATTACAAATGCCGTCGCCATCGACGTCGTTACCGTCGTTCGCGG
CATCGGGGAAGCCGGCACCGGCGGCACAGTCATCGCAACCGTCGCCATCGGCATCGA
>4     
GCGTTCGAAGCAATTGCACGAGACCCAAACAACGAATTGCTGGTTGTTGAACTGGAAAACTCTCTAGTCGGAATGCTTCAAATTACTTATATTCCCTACCTGACACATATTGGCAGTTGGCGTTGTCTTATAGAAGGTGTTCG
AATCCATAGTGACTATCGTGGACGAGGTTTTGGTGAGCAAATGTTCGCACATGCGAT
>5     
GTTTAAGACTAACAGCAATCTGTAAGGACATAGGTGCTGGAGTTGAGGTTAGTCTGGAAGATATGATCTGGGCAGAGAAATTGTCCAAAGCAAACACCGCAGCAAGAGGTATGCTAAACACAGCAAGAAGAATAAGTAATGAT
CCTACTGATTCTTTTCTGAATGAGTTGAATATAGGAGACCCCGACTCAACTCATCAT

输入文件是一个约 8G 的文件，如下所示：

CCAAGATCGCTGGTTGGCGAATCAATTTCATAAACGCCTACGCTTTCAAGGAACGTGTTAAGAATGTTCT
GGCCGAGTTCCTTATGAGACGTTTCGCGTCCCTTAAATCGAATAACGACACGAACCTTGTCGCCGTCATT
AAGAAAACCCTTTGCCTTCTTGGCCTTAATCTGAATATCACGGGTGTCCGTTACAGGTCGCAACTGGATT
TCCTTGACTTCAGAAACAGACTTACGTGAATTCTTCTTGATTTCTTTCTGACGCTTTTCATTTTCATACT
GGAACTTGCCGTAATCAATGATCTTACAAACAGGAATATCACCCTTATCAGAGATCAATACCAAATCAAG
TTCGGCATCAAAAGCGCGATCAAGTGCGTCTTCAATGTCGAGGACCGTTGTTTCTTCACCGTCAACCAAA
CGAATTGTGGAGGACTTGATGTCGTCTCGGGTACTAATTTTATTCACGTATATGTTACTCCTTATGTTGT

任何帮助，将不胜感激。提前致谢。

答案1

这是 awk 中解决方案的实现。该数据是 8GB 的​​伪随机十六进制数字（实际上是大约 12 个手册页的十六进制转换，重复了 3300 次）。大约有 1100 万行，平均每行 725 字节。

这是一个定时执行。

Paul--) ls -l tdHuge.txt
-rw-r--r-- 1 paul paul 8006529300 Dec 24 22:38 tdHuge.txt
Paul--) ./rndSelect
inFile ./tdHuge.txt; Size 8006529300; Count 10000; Lth 200; maxIter 50; Db 1;
Iteration   1 needs  10000
Iteration   2 needs   2712
Overlap   9561: 7663038508 to 7663038657
Iteration   3 needs    728
Iteration   4 needs    195
Iteration   5 needs     50
Iteration   6 needs     11
Iteration   7 needs      2
Required 7 iterations
Reporting 10000 samples

real    2m3.326s
user    0m3.496s
sys 0m10.340s
Paul--) wc Result.txt
  20000   20000 2068894 Result.txt
Paul--) head -n 8 Result.txt | cut -c 1-40
>1
5706C69636174656420696E666F726D6174696F6
>2
20737472696E672028696E207768696368206361
>3
20646F6573206E6F742067657420612068617264
>4
647320616E642073746F7265732E204966207468
Paul--) tail -n 8 Result.txt | cut -c 1-40
>9997
6F7374207369676E69666963616E7420646F7562
>9998
7472696E676F702D73747261746567793D616C67
>9999
865726520736F6D652066696C6573206D7573742
>10000
5726E65642E205768656E20746865202D66206F7
Paul--)

它需要迭代，因为它对文件进行随机探测。如果一个探针与相邻探针或换行符重叠，则该探针将被丢弃，并生成一小批新探针。平均线长度为 725 条线，样本要求为 200 条线，几乎 30% 的探针距离线末端太近而无法接受。我们不知道真实数据的平均行长度——更长的行会提高成功率。

我们也不知道标题行（如 2020 年 12 月 4 日相关问题中所述）是否仍然存在于文件中。但如果每个标题行都小于样本长度 200，则标题行将被丢弃（最好是机缘巧合）。

代码主要是 GNU/awk（最小的 bash），并有一些注释。有很多残留调试，可以通过在选项中设置 Db=0 来隐藏。

#! /bin/bash

#.. Select random non-overlapping character groups from a file.

export LC_ALL="C"

#.. These are the optional values that will need to be edited.
#.. Command options could be used to set these from scripts arguments.

inFile="./tdHuge.txt"
outFile="./Result.txt"
Count=10000     #.. Number of samples.
Lth=200         #.. Length of each sample.
maxIter=50      #.. Prevents excessive attempts.

Size="$( stat -c "%s" "${inFile}" )"
Seed="$( date '+%N' )"
Db=1

#.. Extracts random non-overlapping character groups from a file.

Selector () {

    local Awk='
#.. Seed the random number generation, and show the args being used.
BEGIN {
    NIL = ""; NL = "\n"; SQ = "\047";
    srand (Seed % PROCINFO["pid"]);
    if (Db) printf ("inFile %s; Size %d; Count %d; Lth %d; maxIter %s; Db %s;\n",
        inFile, Size, Count, Lth, maxIter, Db);
    fmtCmd = "dd bs=%d count=1 if=%s iflag=skip_bytes skip=%d status=none";
}
#.. Constructs an array of random file offsets, replacing overlaps.
#.. Existing offsets are indexed from 1 to Count, deleting overlaps.
#.. Additions are indexed from -1 down to -N to avoid clashes.

function Offsets (N, Local, Iter, nSeek, Seek, Text, j) {

    while (N > 0 && Iter < maxIter) {
        ++Iter;
        if (Db) printf ("Iteration %3d needs %6d\n", Iter, N);

        for (j = 1; j <= N; ++j) {
            Seek[-j] = int ((Size - Lth) * rand());
            Text[Seek[-j]] = getSample( Seek[-j], Lth);
            if (Db7) printf ("Added %10d: \n", Seek[-j], Text[Seek[-j]]);
        }
        #.. Reindex in numeric order for overlap checking.
        nSeek = asort (Seek);
        if (Db7) for (j in Seek) printf ("%6d: %10d\n", j, Seek[j]);

        #.. Discard offsets that overlap the next selection.
        N = 0; for (j = 1; j < nSeek; ++j) {
            if (Seek[j] + Lth > Seek[j+1]) {
                if (Db) printf ("Overlap %6d: %10d to %10d\n",
                    j, Seek[j], Seek[j+1]);
                ++N; delete Text[Seek[j]]; delete Seek[j];
            } else if (length (Text[Seek[j]]) < Lth) {
                if (Db7) printf ("Short   %6d: %10d\n",
                    j, Seek[j]);
                ++N; delete Text[Seek[j]]; delete Seek[j];
            }
        }
    }
    if (Iter >= maxIter) {
        printf ("Failed with overlaps after %d iterations\n", Iter);
    } else {
        printf ("Required %d iterations\n", Iter);
        Samples( nSeek, Seek, Text);
    }
}
#.. Returns n bytes from the input file from position p.
function getSample (p, n, Local, cmd, tx) {

    cmd = sprintf (fmtCmd, n, SQ inFile SQ, p);
    if (Db7) printf ("cmd :%s:\n", cmd);
    cmd | getline tx; close (cmd);
    return (tx);
}
#.. Send samples to the output file.
function Samples (nSeek, Seek, Text, Local, j) {

    printf ("Reporting %d samples\n", nSeek);
    for (j = 1; j <= nSeek; ++j) {
        printf (">%d\n%s\n", j, Text[Seek[j]]) > outFile;
    }
    close (outFile);
}
END { Offsets( Count); }
'
    echo | awk -v Size="${Size}" -v inFile="${inFile}" \
        -v outFile="${outFile}" -v Count="${Count}" -v Lth="${Lth}" \
        -v maxIter="${maxIter}" \
        -v Db="${Db}" -v Seed="${Seed}" -f <( printf '%s' "${Awk}" )
}

#.. Test.

    time Selector

答案2

编辑：我将添加另一个答案和实现（可能需要修改）。当前的答案已部分被取代，我可能很快就会删除它。或者也许更新它，作为实际实现的设计说明。

算法的基本框架。它可能不会在统计上发出严格随机的结果，但它似乎确实满足要求。我可能会在 awk 中对此进行原型设计（使用 Cliff RNG，或者可能是从 shuf 读取管道）。

  1. 将 N 设置为要添加到数组 X 的偏移量数量（例如 10000）。

  2. 当 N > 0 时，迭代 2a 到 2d。
    2a.将 N 个随机偏移量块附加到 X（使用现有脚本中的范围）。
    2b.按数字升序对 X 进行排序。
    2c.遍历 X，将每个偏移量与下一个偏移量进行比较，并标记为删除距其后继者 200 个字符以内的任何位置。
    2d.删除标记的偏移量，并将 N 设置为删除的次数。

  3. 此时，我们有 10000 个不重叠的随机偏移（按升序排列）。 （需要测试范围是否足够稀疏，以便第 2 阶段最终退出 - 例如至少 3 * 10000 * 200 个字符。）
    3a。发出一连串10000尾| head 命令（每个偏移量一个）。
    3b.通过 bash 管道传输这些命令。

这需要首先检查范围是否足够稀疏，以便第 2 阶段最终退出——例如至少 3 * 10000 * 200 个字符。可能需要对此进行一些研究。

编辑：检查样本之间的重叠仅需要偏移量，但检查与行尾的重叠需要数据。因此，在 (2a) 处，读取样本本身，并将其存储在按第一个字节的偏移量索引的数组中。检查并删除 (2c) 中的短样本。

这使得整个（3）变得多余，因为我们已经将每个有效样本存储在内存中（并且我们有一个有效偏移量的排序列表来索引它们）。所以我们可以按照偏移量的顺序列出所有样本。

相关内容