DNA 序列状态转移特征提取与 PAC 图绘制

本程序使用 Python 代码实现 DNA 序列的状态转移特征提取，并绘制 PAC 图，用于分析 DNA 序列的特征。程序使用 Biopython 库读取 FASTA 文件，并计算状态转移概率矩阵，然后对矩阵进行标准化并绘制 PAC 图。

代码：

from Bio import SeqIO
import numpy as np
import matplotlib.pyplot as plt
import os

# 读取fasta文件，计算状态转移概率矩阵
def calc_transition_matrix(fasta_file):
    sequences = list(SeqIO.parse(fasta_file, 'fasta'))
    alphabet = ['A', 'C', 'G', 'T']
    states = len(alphabet)
    matrix = np.zeros((states, states))
    for seq in sequences:
        seq_str = str(seq.seq)
        for i in range(len(seq_str) - 1):
            from_state = alphabet.index(seq_str[i])
            to_state = alphabet.index(seq_str[i + 1])
            matrix[from_state, to_state] += 1
    # 将状态转移矩阵转换为状态转移概率矩阵
    matrix = matrix / np.sum(matrix)
    return matrix

# 对矩阵进行z-score标准化
def standardize(matrix):
    standardized_matrix = (matrix - np.mean(matrix, axis=0)) / np.std(matrix, axis=0)
    return standardized_matrix

# 绘制PAC图
def plot_pac(matrix, fig_title):
    cov = np.cov(matrix.T)
    # 对协方差矩阵进行特征值分解，得到特征值和特征向量
    eigenvalues, eigenvectors = np.linalg.eig(cov)
    # 对特征向量进行正交化
    eigenvectors = np.linalg.qr(eigenvectors)[0]
    # 将标准化后的矩阵投影到特征向量上，得到新的矩阵
    new_matrix = np.dot(matrix, eigenvectors)
    # 计算新矩阵的协方差矩阵
    new_cov = np.cov(new_matrix.T)
    # 计算新矩阵的PAC图
    pac = np.zeros_like(new_cov)
    for i in range(new_cov.shape[0]):
        for j in range(new_cov.shape[1]):
            pac[i, j] = new_cov[i, j] / np.sqrt(new_cov[i, i] * new_cov[j, j])
    plt.imshow(pac, cmap='coolwarm')
    plt.colorbar()
    plt.title(fig_title)
    plt.show()

# 测试代码
if __name__ == '__main__':
    # 获取文件夹中的所有fasta文件
    fasta_folder = './FASTA 文件'
    fasta_files = [os.path.join(fasta_folder, f) for f in os.listdir(fasta_folder) if f.endswith('.fasta')]

    # 对每个fasta文件进行状态转移概率矩阵的计算、标准化、PCA分析和PAC图绘制
    for fasta_file in fasta_files:
        # 读取fasta文件，计算状态转移概率矩阵
        matrix = calc_transition_matrix(fasta_file)
        print(matrix)

        # 对矩阵进行标准化
        standardized_matrix = standardize(matrix)
        print(standardized_matrix)

        # 绘制PAC图
        fig_title = os.path.splitext(os.path.basename(fasta_file))[0]
        plot_pac(standardized_matrix, fig_title)

程序使用方法：

将所有 FASTA 文件放在名为 'FASTA 文件' 的文件夹中。
运行程序。
程序将计算每个 FASTA 文件的状态转移概率矩阵，并绘制 PAC 图。

程序特点：

使用 Biopython 库读取 FASTA 文件，方便快捷。
计算状态转移概率矩阵，更准确地反映 DNA 序列的特征。
对矩阵进行标准化，消除不同序列长度的影响。
绘制 PAC 图，直观地显示 DNA 序列的状态转移特征。

程序用途：

分析 DNA 序列的特征。
比较不同 DNA 序列的差异。
识别 DNA 序列的结构和功能。

注意事项：

程序要求安装 Biopython 库。
程序的输入文件必须是 FASTA 格式。
程序生成的 PAC 图可能需要调整颜色和大小才能更好地显示。