如何简化SAS PROC FREQ疾病组合分析并生成热力图？

我来帮你搞定这个问题！用PROC FREQ分析多疾病并发时，动辄生成数千种组合确实让人头大，完全没法人工排查。下面给你一套从筛选有效组合到生成直观热力图的完整方案：

第一步：筛选有意义的疾病并发组合 #

首先我们要把全0、单个疾病的无效组合过滤掉，只保留至少两种疾病并发且有实际病例的组合。可以先把PROC FREQ的结果导出到数据集，再做筛选：

/* 1. 导出PROC FREQ的结果到数据集 */
proc freq data=new; 
    tables HOHT*HOGD*CroD*Psor*Viti*CelD*UlcC*AddD*SluE*Rhea*PerA/out=freq_results(keep=HOHT HOGD CroD Psor Viti CelD UlcC AddD SluE Rhea PerA COUNT) list; 
run;

/* 2. 筛选有效组合：至少2种疾病并发，且病例数>0 */
data filtered_results;
    set freq_results;
    /* 计算当前组合中患病的数量 */
    disease_count = sum(HOHT, HOGD, CroD, Psor, Viti, CelD, UlcC, AddD, SluE, Rhea, PerA);
    /* 过滤掉无效组合 */
    where disease_count >= 2 and COUNT > 0;
    /* 可选：如果想只保留高频组合，再加COUNT阈值，比如COUNT >= 5 */
    /* where disease_count >= 2 and COUNT >= 5; */
run;

这样处理后，剩下的就是真正有价值的并发组合了，数量会大幅减少。

第二步：生成两两疾病的共现矩阵（适合热力图展示） #

如果你的最终目标是热力图，两两疾病的共现关系会比多疾病组合更直观。这里用数组方法快速计算任意两种疾病的并发病例数：

data cooccur_matrix;
    length var1 var2 $30;
    format count 8.;
    set new end=eof;
    /* 把所有疾病变量放进数组 */
    array diseases[*] HOHT HOGD CroD Psor Viti CelD UlcC AddD SluE Rhea PerA;
    /* 初始化计数矩阵，维度和疾病数量一致 */
    retain count_matrix{11,11} 0; 

    /* 遍历所有两两组合，统计并发次数 */
    do i=1 to dim(diseases)-1;
        if diseases[i] = 1 then do j=i+1 to dim(diseases);
            if diseases[j] = 1 then count_matrix[i,j] + 1;
        end;
    end;

    /* 当读取完所有数据后，把矩阵转换成绘图需要的长格式 */
    if eof then do;
        do i=1 to dim(diseases)-1;
            do j=i+1 to dim(diseases);
                var1 = vname(diseases[i]);
                var2 = vname(diseases[j]);
                count = count_matrix[i,j];
                output;
            end;
        end;
    end;
    keep var1 var2 count;
run;

这个数据集里的count就是两种疾病同时发生的病例数，你也可以根据需求转换成并发率（比如除以总病例数，或者其中一种疾病的病例数）。

第三步：用PROC SGPLOT绘制热力图 #

现在用PROC SGPLOT的HEATMAPCATEGORICAL语句生成热力图，还可以加上数值标签让结果更清晰：

proc sgplot data=cooccur_matrix;
    /* 绘制热力图，颜色随并发数量深浅变化 */
    heatmapcategorical x=var1 y=var2 colorresponse=count / 
        colormodel=(white lightcoral darkred) /* 自定义颜色渐变 */
        display=(fill outline)
        name='heat';
    /* 在每个单元格添加并发数量的标签 */
    text x=var1 y=var2 text=count / textattrs=(size=10 color=white) position=center;
    title "疾病并发情况热力图";
    /* 隐藏坐标轴标签，让图表更简洁 */
    xaxis display=(nolabel) fitpolicy=rotate;
    yaxis display=(nolabel);
run;

如果想要对称的热力图（比如A-B和B-A显示相同的数值），可以再补充对称的行到cooccur_matrix数据集中，或者调整绘图参数。

内容的提问来源于stack exchange，提问作者ybao