如何使用Python正则从fastq文件中提取指定精确核苷酸序列

• 未遵循fastq的4行结构规则：fastq每4行对应一条测序读长，仅第2行是核苷酸序列，原代码逐行扫描所有行，无法准确定位到序列行，还会错误处理质量值行、头行等非序列内容。
• 正则方法使用错误：re.match()仅从字符串开头位置匹配，目标序列可能出现在序列行的任意位置，应该用re.search()；另外待匹配的是固定序列，不需要调用正则，直接用in判断效率更高。
• 未做输入预处理：测试fastq每行开头有大量空格，需要先去除首尾空白再做判断，否则头行匹配、序列匹配都会失败。
• 输出内容错误：匹配成功后原代码直接打印re.match()返回的匹配对象，无法输出可读的序列或读长信息。

注：你当前提供的测试文件中所有序列都不包含目标核苷酸序列，替换为带目标序列的测试文件即可得到输出。

target_seq = "GATCGGAAGAGCTCGTATGCCGTCTTCTGCTTGAAA"
line_counter = 0
current_header = ""

with open('last_mock.fastq', 'r', encoding='utf-8') as rf:
    for line in rf:
        cleaned_line = line.strip()
        mod = line_counter % 4
        if mod == 0:
            # 存储当前读长的头行信息
            current_header = cleaned_line
            is_valid_read = cleaned_line.startswith('@')
        elif mod == 1 and is_valid_read:
            # 仅对序列行做匹配判断
            if target_seq in cleaned_line:
                # 可按需调整输出，这里同时输出头行和序列
                print(current_header)
                print(cleaned_line)
        line_counter += 1

• 如果需要提取目标序列在序列行中的位置，可以用cleaned_line.find(target_seq)获取起始索引。
• 如果要严格校验fastq格式，可以增加对第3行+开头的判断，避免异常格式文件的干扰。

内容的提问来源于stack exchange，提问作者kira.99