Pandas中value_counts、.loc、.iloc、np.linspace()用法解释

根据这个代码画一个流程图import numpy as np import pandas as pd from sklearn.cluster import KMeans from sklearn.preprocessing import StandardScaler import matplotlib.pyplot as plt import seaborn as sns import os import sys # 设置中文字体 plt.rcParams["font.family"] = "SimHei" plt.rcParams["axes.unicode_minus"] = False # 目标目录 TARGET_DIR = "/mnt" # 创建目标目录 def create_target_directory(): if not os.path.exists(TARGET_DIR): try: os.makedirs(TARGET_DIR) print(f"目录 {TARGET_DIR} 创建成功") except OSError as e: print(f"创建目录失败: {e}") sys.exit(1) # 数据读取与预处理 def load_and_preprocess_data(file_path): """加载并预处理航空数据""" try: # 尝试使用utf-8编码读取 airline_data = pd.read_csv(file_path, encoding="utf-8") print("使用utf-8编码成功读取文件") except UnicodeDecodeError: try: # 尝试使用gb18030编码读取 airline_data = pd.read_csv(file_path, encoding="gb18030") print("使用gb18030编码成功读取文件") except UnicodeDecodeError: # 使用更通用的错误处理方式 airline_data = pd.read_csv(file_path, encoding="utf-8", errors="replace") print("使用错误替换模式读取文件，可能存在部分字符替换") print('原始数据的形状为：', airline_data.shape) # 去除票价为空的记录 exp1 = airline_data["SUM_YR_1"].notnull() exp2 = airline_data["SUM_YR_2"].notnull() exp = exp1 & exp2 airline_notnull = airline_data.loc[exp, :] print('删除缺失记录后数据的形状为：', airline_notnull.shape) # 只保留票价非零的，或者平均折扣率不为0且总飞行公里数大于0的记录 index1 = airline_notnull['SUM_YR_1'] != 0 index2 = airline_notnull['SUM_YR_2'] != 0 index3 = (airline_notnull['SEG_KM_SUM'] > 0) & (airline_notnull['avg_discount'] != 0) airline = airline_notnull[(index1 | index2) & index3] print('删除异常记录后数据的形状为：', airline.shape) # 检查并清洗FFP_TIER和AGE字段 airline = airline[airline['FFP_TIER'].notnull()] airline = airline[pd.to_numeric(airline['FFP_TIER'], errors='coerce').notnull()] airline['FFP_TIER'] = pd.to_numeric(airline['FFP_TIER']) airline = airline[airline['AGE'].notnull()] airline = airline[pd.to_numeric(airline['AGE'], errors='coerce').notnull()] airline['AGE'] = pd.to_numeric(airline['AGE']) print('删除异常记录后数据的形状为：', airline.shape) return airline # 特征工程 def feature_engineering(airline_data): """构建LRFMC特征""" # 选取需求特征 airline_selection = airline_data[["FFP_DATE", "LOAD_TIME", "FLIGHT_COUNT", "LAST_TO_END", "avg_discount", "SEG_KM_SUM"]] # 构建L特征 (会员关系长度，单位：月) print("正在计算会员关系长度(L)...") L = pd.to_datetime(airline_selection["LOAD_TIME"]) - pd.to_datetime(airline_selection["FFP_DATE"]) L = L.astype("str").str.split().str[0] L = L.astype("int") / 30 # 合并特征并命名为LRFMC print("正在合并特征...") airline_features = pd.concat([L, airline_selection.iloc[:, 2:]], axis=1) airline_features.columns = ['L', 'R', 'F', 'M', 'C'] # 重命名列名 print('构建的LRFMC特征前5行为：\n', airline_features.head()) return airline_features # 数据标准化 def standardize_features(features): """标准化特征数据并保存到目标目录""" print("正在标准化特征数据...") scaler = StandardScaler() data = scaler.fit_transform(features) np.savez(os.path.join(TARGET_DIR, "airline_scale.npz"), data) print('标准化后LRFMC五个特征前5行为：\n', data[:5, :]) return data, scaler # KMeans聚类 def kmeans_clustering(data, n_clusters=5, random_state=123): """使用KMeans进行客户分群""" print(f"正在进行KMeans聚类，聚类中心数: {n_clusters}") # 构建模型 kmeans_model = KMeans(n_clusters=n_clusters, init='k-means++', n_init='auto', random_state=random_state) fit_kmeans = kmeans_model.fit(data) # 模型训练 # 查看聚类结果 cluster_centers = kmeans_model.cluster_centers_ labels = kmeans_model.labels_ # 统计不同类别样本的数目 cluster_counts = pd.Series(labels).value_counts().sort_index() print('最终每个类别的数目为：\n', cluster_counts) return kmeans_model, cluster_centers, labels, cluster_counts # 计算最佳聚类数 def calculate_optimal_clusters(data, max_clusters=10): """使用肘部法则确定最佳聚类数并保存图片""" print("正在计算最佳聚类数...") wss = [] for k in range(1, max_clusters + 1): kmeans = KMeans(n_clusters=k, init='k-means++', n_init='auto', random_state=123) kmeans.fit(data) wss.append(kmeans.inertia_) # 可视化肘部曲线并保存到目标目录 plt.figure(figsize=(10, 6)) plt.plot(range(1, max_clusters + 1), wss, 'bo-') plt.xlabel('聚类中心数(k)') plt.ylabel('簇内误差平方和(WSS)') plt.title('肘部法则确定最佳聚类数') plt.savefig(os.path.join(TARGET_DIR, "elbow_method.png")) plt.close() print("肘部曲线已保存至", os.path.join(TARGET_DIR, "elbow_method.png")) return wss # 可视化聚类结果 def visualize_clusters(cluster_centers, cluster_counts, feature_names=['L', 'R', 'F', 'M', 'C']): """可视化聚类中心和分布并保存图片""" n_clusters = cluster_centers.shape[0] n_features = cluster_centers.shape[1] # 1. 雷达图保存到目标目录 angles = np.linspace(0, 2 * np.pi, n_features, endpoint=False).tolist() angles += angles[:1] # 闭合 fig, ax = plt.subplots(figsize=(10, 10), subplot_kw=dict(polar=True)) for i in range(n_clusters): values = cluster_centers[i].tolist() values += values[:1] # 闭合 ax.plot(angles, values, 'o-', linewidth=2, label=f'Cluster {i} ({cluster_counts[i]}客户)') ax.fill(angles, values, alpha=0.25) ax.set_thetagrids(np.degrees(angles[:-1]), feature_names) ax.set_ylim(cluster_centers.min() - 0.5, cluster_centers.max() + 0.5) plt.legend(loc='upper right') plt.title('客户聚类中心雷达图') plt.savefig(os.path.join(TARGET_DIR, "customer_clusters_radar.png")) plt.close() # 2. 柱状图保存到目标目录 plt.figure(figsize=(10, 6)) bars = plt.bar(range(n_clusters), cluster_counts.values, tick_label=[f'Cluster {i}' for i in range(n_clusters)]) plt.xlabel('客户群') plt.ylabel('客户数量') plt.title('客户群分布') # 添加数值标签 for bar in bars: height = bar.get_height() plt.text(bar.get_x() + bar.get_width() / 2., height, f'{height}', ha='center', va='bottom') plt.savefig(os.path.join(TARGET_DIR, "customer_clusters_distribution.png")) plt.close() print("聚类可视化结果已保存至", TARGET_DIR) # 分析并解释聚类结果（新增营销策略建议） def analyze_clusters(features, labels, scaler, n_clusters=5): """分析并解释聚类结果，添加营销策略建议""" # 还原聚类中心到原始特征空间 cluster_centers_original = scaler.inverse_transform(np.array([ features.iloc[labels == i].mean().values for i in range(n_clusters) ])) # 创建聚类中心DataFrame center_df = pd.DataFrame(cluster_centers_original, columns=['L', 'R', 'F', 'M', 'C']) center_df.index.name = 'Cluster' # 计算各聚类占比 cluster_counts = pd.Series(labels).value_counts().sort_index() total = len(labels) center_df['客户数量'] = cluster_counts center_df['占比(%)'] = (cluster_counts / total * 100).round(2) # 添加营销策略建议 high_value_threshold = { 'M': features['M'].mean() * 1.5, # 飞行里程高于均值1.5倍 'F': features['F'].mean() * 1.2 # 飞行频率高于均值1.2倍 } active_threshold = 30 # 最近消费间隔小于30天视为活跃 # 生成营销策略 center_df['营销策略'] = '' for i in range(n_clusters): m_value = center_df.loc[i, 'M'] f_value = center_df.loc[i, 'F'] r_value = center_df.loc[i, 'R'] if m_value > high_value_threshold['M'] and f_value > high_value_threshold['F']: # 高价值高忠诚客户 center_df.loc[i, '营销策略'] = """ ★ 顶级服务策略 ★ - 专属客户经理一对一服务 - 优先升舱/贵宾厅权益 - 高价值积分兑换（如免费国际航班） - 生日/节日专属礼品 - 邀请参与会员体验活动 """ elif r_value < active_threshold and m_value > features['M'].mean(): # 活跃中价值客户 center_df.loc[i, '营销策略'] = """ ★ 活跃激励策略 ★ - 推送个性化套餐（如商务舱折扣） - 飞行里程加倍活动 - 邀请参与会员专属活动 - 累计消费达标赠礼 - 推荐联名信用卡 """ elif m_value < features['M'].mean() and f_value < features['F'].mean(): # 低价值潜在客户 center_df.loc[i, '营销策略'] = """ ★ 唤醒激活策略 ★ - 发放新客/复购优惠券（如满减机票） - 推荐热门旅游航线套餐 - 首次飞行赠送基础积分 - 社交分享获额外折扣 - 发送限时特价提醒 """ else: # 其他类型客户（如低频高折扣客户） center_df.loc[i, '营销策略'] = """ ★ 定向培育策略 ★ - 发送定制化服务通知（如低价预警） - 团体票/家庭票优惠推荐 - 积分过期提醒与补签机制 - 调研客户需求优化服务 - 推荐会员等级晋升计划 """ # 保存到目标目录 center_df.to_csv(os.path.join(TARGET_DIR, "cluster_analysis.csv"), index=True) # 打印聚类分析结果 print("\n===== 客户分群分析结果 =====") print(center_df[['客户数量', '占比(%)', 'L', 'R', 'F', 'M', 'C', '营销策略']]) # 客户群特征描述 print("\n===== 客户群特征核心指标 =====") for i in range(n_clusters): print(f"\n客户群 {i} - 核心特征：") print( f" ● 里程(M)：{center_df.loc[i, 'M']:.1f}公里 ({'高于' if center_df.loc[i, 'M'] > features['M'].mean() else '低于'}均值)") print( f" ● 频率(F)：{center_df.loc[i, 'F']:.1f}次 ({'活跃' if center_df.loc[i, 'F'] > features['F'].mean() else '低频'})") print( f" ● 最近消费(R)：{center_df.loc[i, 'R']:.1f}天 ({'近期活跃' if center_df.loc[i, 'R'] < active_threshold else '长期未活跃'})") return center_df # 绘制性别比例饼图 def plot_gender_pie(airline_data): gender_counts = airline_data['GENDER'].value_counts() plt.figure(figsize=(6, 6)) plt.pie(gender_counts, labels=gender_counts.index, autopct='%1.1f%%', startangle=90) plt.title('性别比例') plt.savefig(os.path.join(TARGET_DIR, "gender_pie_chart.png")) plt.close() print("性别比例饼图已保存至", os.path.join(TARGET_DIR, "gender_pie_chart.png")) # 绘制热力图 def plot_heatmap(airline_data): # 选取关键列 key_columns = ['FLIGHT_COUNT', 'LAST_TO_END', 'SEG_KM_SUM', 'EXCHANGE_COUNT', 'Points_Sum'] numeric_data = airline_data[key_columns] corr_matrix = numeric_data.corr() plt.figure(figsize=(10, 8)) sns.heatmap(corr_matrix, annot=True, cmap='coolwarm', vmin=-1, vmax=1) plt.title('数值型变量相关性热力图') plt.savefig(os.path.join(TARGET_DIR, "heatmap.png")) plt.close() print("热力图已保存至", os.path.join(TARGET_DIR, "heatmap.png")) # 绘制客户年龄与累计飞行里程关系的散点图 def plot_age_mileage_scatter(airline_data): plt.figure(figsize=(10, 6)) plt.scatter(airline_data['AGE'], airline_data['SEG_KM_SUM']) plt.xlabel('客户年龄') plt.xticks(rotation=45) plt.ylabel('累计飞行里程') plt.title('客户年龄与累计飞行里程关系散点图') plt.savefig(os.path.join(TARGET_DIR, "age_mileage_scatter.png")) plt.close() print("客户年龄与累计飞行里程关系散点图已保存至", os.path.join(TARGET_DIR, "age_mileage_scatter.png")) def plot_member_level_age(airline_data): if 'FFP_TIER' not in airline_data.columns or 'AGE' not in airline_data.columns: print("错误：数据中缺少会员等级(FFP_TIER)或年龄(AGE)字段！") return # 定义与实际数据匹配的映射关系 tier_mapping = { 4: '银卡会员', 5: '金卡会员', 6: '白金卡会员', } airline_data['会员等级'] = airline_data['FFP_TIER'].map(tier_mapping).fillna('未知等级') # 计算各等级的平均年龄和人数 age_avg = airline_data.groupby('会员等级')['AGE'].agg(['mean', 'count']).reset_index() age_avg.columns = ['会员等级', '平均年龄', '人数'] # 按等级顺序排序 order = ['银卡会员', '金卡会员', '白金卡会员', '未知等级'] age_avg['会员等级'] = pd.Categorical(age_avg['会员等级'], categories=order, ordered=True) age_avg = age_avg.sort_values('会员等级') # 绘制柱状图 plt.figure(figsize=(10, 6)) bars = plt.bar(age_avg['会员等级'], age_avg['平均年龄'], width=0.6, color=['#4CAF50', '#2196F3', '#FF9800', '#9E9E9E']) plt.xlabel('会员等级') plt.ylabel('平均年龄（岁）') plt.title('不同会员等级（FFP_TIER）平均年龄分布') # 添加数值标签（直接使用索引获取数据） for i, bar in enumerate(bars): height = bar.get_height() # 通过索引直接获取对应等级的人数 count = age_avg.iloc[i]['人数'] plt.text(bar.get_x() + bar.get_width()/2, height, f'平均：{height:.1f}岁\n({count}人)', ha='center', va='bottom', fontsize=10) plt.grid(axis='y', linestyle='--', alpha=0.7) plt.xticks(rotation=0) plt.tight_layout() plt.savefig(os.path.join(TARGET_DIR, "member_level_age.png")) plt.close() print("会员等级平均年龄柱状图已保存至", os.path.join(TARGET_DIR, "member_level_age.png")) # 主函数 def main(file_path="air_data.csv", n_clusters=5): """主函数：执行完整的客户分群分析流程""" print(f"===== 航空客户分群分析系统 =====") # 创建目标目录 create_target_directory() # 数据加载与预处理 print("\n===== 数据加载与预处理 =====") airline_data = load_and_preprocess_data(file_path) # 特征工程 print("\n===== 特征工程 =====") features = feature_engineering(airline_data) # 数据标准化 print("\n===== 数据标准化 =====") standardized_data, scaler = standardize_features(features) # 计算最佳聚类数 calculate_optimal_clusters(standardized_data, max_clusters=10) # KMeans聚类 print(f"\n===== KMeans聚类 (k={n_clusters}) =====") kmeans_model, centers, labels, counts = kmeans_clustering( standardized_data, n_clusters=n_clusters ) # 可视化聚类结果 print("\n===== 可视化分析 =====") visualize_clusters(centers, counts) # 绘制性别比例饼图 print("\n===== 绘制性别比例饼图 =====") plot_gender_pie(airline_data) # 绘制热力图 print("\n===== 绘制热力图 =====") plot_heatmap(airline_data) # 绘制客户年龄与累计飞行里程关系的散点图 print("\n===== 绘制客户年龄与累计飞行里程关系的散点图 =====") plot_age_mileage_scatter(airline_data) # 调用绘制不同会员等级的平均年龄柱状图函数 print("\n===== 绘制会员等级平均年龄分布 =====") plot_member_level_age(airline_data) # 分析并解释聚类结果 print("\n===== 结果分析 =====") cluster_analysis = analyze_clusters(features, labels, scaler, n_clusters) # 保存完整聚类结果到目标目录 result_df = pd.concat([ features.reset_index(drop=True), pd.Series(labels, name='Cluster') ], axis=1) result_df.to_csv(os.path.join(TARGET_DIR, "customer_clusters.csv"), index=False) print(f"\n===== 分析完成 =====") print(f"所有结果已保存至 {TARGET_DIR} 目录：") print("- 标准化数据: airline_scale.npz") print("- 聚类结果: customer_clusters.csv") print("- 聚类分析: cluster_analysis.csv") print("- 可视化图表: customer_clusters_radar.png, customer_clusters_distribution.png, elbow_method.png, gender_pie_chart.png, heatmap.png") print("- 客户年龄与累计飞行里程: age_mileage_scatter.png") print("- 会员等级平均年龄: member_level_age.png") if __name__ == "__main__": main()

import pandas as pd import numpy as np import matplotlib.pyplot as plt import seaborn as sns from scipy import stats import os from statsmodels.stats.multicomp import pairwise_tukeyhsd from statsmodels.formula.api import ols import statsmodels.api as sm from scipy.stats import chi2_contingency, fisher_exact import warnings warnings.filterwarnings('ignore') # 设置matplotlib参数，使图表更美观，支持中文显示 plt.rcParams['font.sans-serif'] = ['SimHei'] # 设置中文字体 plt.rcParams['axes.unicode_minus'] = False # 解决保存图像时负号'-'显示为方块的问题 plt.rcParams['figure.dpi'] = 300 # 设置图像DPI为300 plt.rcParams['savefig.dpi'] = 300 # 设置保存图像的DPI为300 # 创建保存结果的文件夹 if not os.path.exists('问题1结果'): os.makedirs('问题1结果') # 读取数据 file_path = '附件1：两个医院临床受试者及抑郁症的基本数据.xlsx' # 读取两个子表 df1 = pd.read_excel(file_path, sheet_name='一院临床受试者及抑郁症的基本数据', header=1) df2 = pd.read_excel(file_path, sheet_name='二院临床受试者及抑郁症的基本数据', header=1) # 重命名列以匹配实际数据 columns = ['序号', '组别', '年龄_岁', '未婚', '已婚', '离异', '丧偶', '无', '使用过抗抑郁药', '其它', '轻度', '中度', '重度'] df1.columns = columns df2.columns = columns # 添加医院标识 df1['医院'] = '医院1' df2['医院'] = '医院2' # 合并数据集 combined_df = pd.concat([df1, df2], ignore_index=True) # 数据预处理 # 将婚姻状况转换为单列分类变量 combined_df['婚姻状况'] = combined_df[['未婚', '已婚', '离异', '丧偶']].idxmax(axis=1) # 将既往抗抑郁药使用情况转换为单列分类变量 combined_df['既往用药'] = combined_df[['无', '使用过抗抑郁药', '其它']].idxmax(axis=1) # 处理抑郁程度列 for col in ['轻度', '中度', '重度']: combined_df[col] = combined_df[col].astype(str) combined_df[col] = combined_df[col].apply(lambda x: 1 if x == '1' else 0) combined_df['抑郁程度'] = combined_df[['轻度', '中度', '重度']].idxmax(axis=1) # 转换组别为字符串类型 combined_df['组别'] = combined_df['组别'].astype(str) # 添加药物名称对应关系 drug_map = {'1': '药物A', '2': '药物B', '3': '药物C'} combined_df['药物名称'] = combined_df['组别'].map(drug_map) # 年龄分析 # 描述性统计 age_stats = combined_df.groupby('组别')['年龄_岁'].describe() age_stats.to_excel('问题1结果/问题1_年龄描述性统计.xlsx') # 绘制年龄分布箱线图 plt.figure(figsize=(10, 6)) sns.boxplot(x='药物名称', y='年龄_岁', data=combined_df, palette='Set3') plt.title('各组受试者年龄分布箱线图') plt.xlabel('药物名称') plt.ylabel('年龄（岁）') plt.grid(True, linestyle='--', alpha=0.7) plt.tight_layout() plt.savefig('问题1结果/问题1_年龄分布箱线图.png') plt.show() # 绘制年龄分布核密度图 plt.figure(figsize=(10, 6)) for group in sorted(combined_df['药物名称'].unique()): group_data = combined_df[combined_df['药物名称'] == group]['年龄_岁'] sns.kdeplot(group_data, label=group, shade=True, alpha=0.5) plt.title('各组受试者年龄分布核密度图') plt.xlabel('年龄（岁）') plt.ylabel('密度') plt.legend() plt.grid(True, linestyle='--', alpha=0.7) plt.tight_layout() plt.savefig('问题1结果/问题1_年龄分布核密度图.png') plt.show() # ANOVA分析检验年龄差异 formula = '年龄_岁 ~ C(组别)' model = ols(formula, data=combined_df).fit() anova_table = sm.stats.anova_lm(model, typ=2) print("年龄ANOVA分析结果：") print(anova_table) # 计算效应量 eta_squared = anova_table['sum_sq'][0] / (anova_table['sum_sq'][0] + anova_table['sum_sq'][1]) print(f"年龄方差分析效应量 (Eta Squared): {eta_squared:.4f}") # 事后检验 tukey_results = pairwise_tukeyhsd(combined_df['年龄_岁'], combined_df['组别']) print("年龄Tukey HSD事后检验结果:") print(tukey_results) # 将Tukey结果转换为DataFrame并保存 tukey_df = pd.DataFrame(data=tukey_results._results_table.data[1:], columns=tukey_results._results_table.data[0]) tukey_df.to_excel('问题1结果/问题1_年龄Tukey事后检验.xlsx', index=False) # 婚姻状况分析 marriage_counts = pd.crosstab(combined_df['药物名称'], combined_df['婚姻状况']) marriage_percentages = pd.crosstab(combined_df['药物名称'], combined_df['婚姻状况'], normalize='index') * 100 marriage_counts.to_excel('问题1结果/问题1_婚姻状况分布频数.xlsx') marriage_percentages.to_excel('问题1结果/问题1_婚姻状况分布百分比.xlsx') # 绘制婚姻状况分布堆积柱状图 plt.figure(figsize=(10, 6)) marriage_counts.plot(kind='bar', stacked=True, colormap='Paired') plt.title('各药物组受试者婚姻状况分布') plt.xlabel('药物名称') plt.ylabel('人数') plt.legend(title='婚姻状况') plt.grid(axis='y', linestyle='--', alpha=0.7) plt.tight_layout() plt.savefig('问题1结果/问题1_婚姻状况分布堆积柱状图.png') plt.show() # 绘制婚姻状况分布热图 plt.figure(figsize=(10, 6)) sns.heatmap(marriage_percentages, annot=True, fmt='.1f', cmap='YlGnBu') plt.title('各药物组受试者婚姻状况分布百分比热图') plt.xlabel('婚姻状况') plt.ylabel('药物名称') plt.tight_layout() plt.savefig('问题1结果/问题1_婚姻状况分布热图.png') plt.show() # 绘制婚姻状况分布堆积柱状图（带比例标注） fig, ax = plt.subplots(figsize=(10, 6)) marriage_counts.plot(kind='bar', stacked=True, colormap='Paired', ax=ax) for idx, bar in enumerate(ax.patches): if idx % 4 == 0: height = bar.get_height() width = bar.get_width() x = bar.get_x() + width / 2 y = height + 0.1 ax.text(x, y, f'{height/sum(marriage_counts.iloc[int(idx/4)])/4:.1%}', ha='center', fontsize=9) plt.title('各药物组受试者婚姻状况分布（带比例标注）') plt.xlabel('药物名称') plt.ylabel('人数') plt.legend(title='婚姻状况') plt.grid(axis='y', linestyle='--', alpha=0.7) plt.tight_layout() plt.savefig('问题1结果/问题1_婚姻状况分布堆积柱状图优化.png') plt.show() # 卡方检验婚姻状况差异 chi2, p, dof, expected = chi2_contingency(marriage_counts) print("婚姻状况卡方检验结果：") print(f"卡方值: {chi2:.4f}, p值: {p:.4f}, 自由度: {dof}") # 计算Cramer's V效应量 n = marriage_counts.sum().sum() cramer_v = np.sqrt(chi2 / (n * (min(marriage_counts.shape) - 1))) print(f"Cramer's V效应量: {cramer_v:.4f}") # 保存婚姻状况卡方检验结果 chi2_results = pd.DataFrame({ '统计量': ['卡方值', 'p值', '自由度', 'Cramer\'s V'], '数值': [chi2, p, dof, cramer_v] }) chi2_results.to_excel('问题1结果/问题1_婚姻状况卡方检验结果.xlsx', index=False) # 既往抗抑郁药使用情况分析 med_history_counts = pd.crosstab(combined_df['药物名称'], combined_df['既往用药']) med_history_percentages = pd.crosstab(combined_df['药物名称'], combined_df['既往用药'], normalize='index') * 100 med_history_counts.to_excel('问题1结果/问题1_既往用药分布频数.xlsx') med_history_percentages.to_excel('问题1结果/问题1_既往用药分布百分比.xlsx') # 绘制既往用药分布堆积柱状图 plt.figure(figsize=(10, 6)) med_history_counts.plot(kind='bar', stacked=True, colormap='Set2') plt.title('各药物组受试者既往抗抑郁药使用情况分布') plt.xlabel('药物名称') plt.ylabel('人数') plt.legend(title='既往用药情况') plt.grid(axis='y', linestyle='--', alpha=0.7) plt.tight_layout() plt.savefig('问题1结果/问题1_既往用药分布堆积柱状图.png') plt.show() # 绘制既往用药分布饼图组合 fig, axes = plt.subplots(1, 3, figsize=(15, 5)) groups = sorted(combined_df['药物名称'].unique()) for i, group in enumerate(groups): group_data = med_history_counts.loc[group] axes[i].pie(group_data, labels=group_data.index, autopct='%1.1f%%', startangle=90, colors=sns.color_palette('Set2', n_colors=len(group_data))) axes[i].set_title(f'{group} 既往用药情况') plt.tight_layout() plt.savefig('问题1结果/问题1_既往用药分布饼图组合.png') plt.show() # 绘制既往用药分布堆积柱状图（带比例标注） fig, ax = plt.subplots(figsize=(10, 6)) med_history_counts.plot(kind='bar', stacked=True, colormap='Set2', ax=ax) for idx, bar in enumerate(ax.patches): if idx % 3 == 0: height = bar.get_height() width = bar.get_width() x = bar.get_x() + width / 2 y = height + 0.1 ax.text(x, y, f'{height/sum(med_history_counts.iloc[int(idx/3)])/3:.1%}', ha='center', fontsize=9) plt.title('各药物组受试者既往抗抑郁药使用情况分布（带比例标注）') plt.xlabel('药物名称') plt.ylabel('人数') plt.legend(title='既往用药情况') plt.grid(axis='y', linestyle='--', alpha=0.7) plt.tight_layout() plt.savefig('问题1结果/问题1_既往用药分布堆积柱状图优化.png') plt.show() # 卡方检验既往用药差异 chi2, p, dof, expected = chi2_contingency(med_history_counts) print("既往抗抑郁药使用情况卡方检验结果：") print(f"卡方值: {chi2:.4f}, p值: {p:.4f}, 自由度: {dof}") # 计算Cramer's V效应量 n = med_history_counts.sum().sum() cramer_v = np.sqrt(chi2 / (n * (min(med_history_counts.shape) - 1))) print(f"Cramer's V效应量: {cramer_v:.4f}") # 保存既往用药卡方检验结果 chi2_results = pd.DataFrame({ '统计量': ['卡方值', 'p值', '自由度', 'Cramer\'s V'], '数值': [chi2, p, dof, cramer_v] }) chi2_results.to_excel('问题1结果/问题1_既往用药卡方检验结果.xlsx', index=False) # 初始抑郁程度分析 depression_counts = pd.crosstab(combined_df['药物名称'], combined_df['抑郁程度']) depression_percentages = pd.crosstab(combined_df['药物名称'], combined_df['抑郁程度'], normalize='index') * 100 depression_counts.to_excel('问题1结果/问题1_抑郁程度分布频数.xlsx') depression_percentages.to_excel('问题1结果/问题1_抑郁程度分布百分比.xlsx') # 绘制抑郁程度分布堆积柱状图 plt.figure(figsize=(10, 6)) depression_counts.plot(kind='bar', stacked=True, colormap='viridis') plt.title('各药物组受试者初始抑郁程度分布') plt.xlabel('药物名称') plt.ylabel('人数') plt.legend(title='抑郁程度') plt.grid(axis='y', linestyle='--', alpha=0.7) plt.tight_layout() plt.savefig('问题1结果/问题1_抑郁程度分布堆积柱状图.png') plt.show() # 绘制抑郁程度分布马赛克图 plt.figure(figsize=(10, 6)) from statsmodels.graphics.mosaicplot import mosaic depression_counts_df = pd.DataFrame(depression_counts) depression_index = depression_counts_df.index.tolist() depression_columns = depression_counts_df.columns.tolist() mosaic_data = {} for idx in depression_index: for col in depression_columns: mosaic_data[(idx, col)] = depression_counts.loc[idx, col] mosaic(mosaic_data, gap=0.02, title='各药物组受试者初始抑郁程度分布马赛克图', properties=lambda key: {'color': plt.cm.viridis(depression_columns.index(key[1]) / len(depression_columns))}) plt.tight_layout() plt.savefig('问题1结果/问题1_抑郁程度分布马赛克图.png') plt.show() # 绘制抑郁程度分布堆积柱状图（带比例标注） fig, ax = plt.subplots(figsize=(10, 6)) depression_counts.plot(kind='bar', stacked=True, colormap='viridis', ax=ax) for idx, bar in enumerate(ax.patches): if idx % 3 == 0: height = bar.get_height() width = bar.get_width() x = bar.get_x() + width / 2 y = height + 0.1 ax.text(x, y, f'{height/sum(depression_counts.iloc[int(idx/3)])/3:.1%}', ha='center', fontsize=9) plt.title('各药物组受试者初始抑郁程度分布（带比例标注）') plt.xlabel('药物名称') plt.ylabel('人数') plt.legend(title='抑郁程度') plt.grid(axis='y', linestyle='--', alpha=0.7) plt.tight_layout() plt.savefig('问题1结果/问题1_抑郁程度分布堆积柱状图优化.png') plt.show() # 卡方检验抑郁程度差异 chi2, p, dof, expected = chi2_contingency(depression_counts) print("初始抑郁程度卡方检验结果：") print(f"卡方值: {chi2:.4f}, p值: {p:.4f}, 自由度: {dof}") # 计算Cramer's V效应量 n = depression_counts.sum().sum() cramer_v = np.sqrt(chi2 / (n * (min(depression_counts.shape) - 1))) print(f"Cramer's V效应量: {cramer_v:.4f}") # 保存抑郁程度卡方检验结果 chi2_results = pd.DataFrame({ '统计量': ['卡方值', 'p值', '自由度', 'Cramer\'s V'], '数值': [chi2, p, dof, cramer_v] }) chi2_results.to_excel('问题1结果/问题1_抑郁程度卡方检验结果.xlsx', index=False) # 综合基线特征雷达图 categories = ['平均年龄', '未婚比例(%)', '无既往用药比例(%)', '重度抑郁比例(%)'] drug_labels = sorted(combined_df['药物名称'].unique()) # 计算各组的值 baseline_values = {} for drug in drug_labels: group_data = combined_df[combined_df['药物名称'] == drug] age_mean = group_data['年龄_岁'].mean() unmarried_pct = marriage_percentages.loc[drug, '未婚'] if drug in marriage_percentages.index else 0 no_med_history_pct = med_history_percentages.loc[drug, '无'] if drug in med_history_percentages.index else 0 severe_depression_pct = depression_percentages.loc[drug, '重度'] if drug in depression_percentages.index else 0 baseline_values[drug] = [age_mean, unmarried_pct, no_med_history_pct, severe_depression_pct] # 绘制雷达图 angles = np.linspace(0, 2 * np.pi, len(categories), endpoint=False).tolist() angles += angles[:1] # 闭合图形 fig, ax = plt.subplots(figsize=(10, 8), subplot_kw=dict(polar=True)) for drug, values in baseline_values.items(): values += values[:1] # 闭合图形 ax.plot(angles, values, linewidth=2, label=f'{drug}') ax.fill(angles, values, alpha=0.1) # 修正变量名 ax.set_thetagrids(np.degrees(angles[:-1]), categories) ax.set_ylim(0, 100) ax.set_title('各药物组基线特征雷达图比较', size=15) ax.legend(loc='upper right', bbox_to_anchor=(0.1, 0.1)) plt.tight_layout() plt.savefig('问题1结果/问题1_基线特征雷达图.png') plt.show() # 基线特征总结表格 age_summary = combined_df.groupby('组别')['年龄_岁'].agg(['mean', 'std']).round(2) age_summary.columns = ['平均年龄', '年龄标准差'] marriage_summary = marriage_percentages.round(1) med_history_summary = med_history_percentages.round(1) depression_summary = depression_percentages.round(1) summary_columns = pd.MultiIndex.from_tuples([ ('年龄', '平均值'), ('年龄', '标准差'), ('婚姻状况', '未婚%'), ('婚姻状况', '已婚%'), ('婚姻状况', '离异%'), ('婚姻状况', '丧偶%'), ('既往用药', '无%'), ('既往用药', '使用过%'), ('既往用药', '其它%'), ('抑郁程度', '轻度%'), ('抑郁程度', '中度%'), ('抑郁程度', '重度%') ]) baseline_summary = pd.DataFrame(index=sorted(combined_df['组别'].unique()), columns=summary_columns) for group in sorted(combined_df['组别'].unique()): baseline_summary.loc[group, ('年龄', '平均值')] = age_summary.loc[group, '平均年龄'] baseline_summary.loc[group, ('年龄', '标准差')] = age_summary.loc[group, '年龄标准差'] for status in marriage_percentages.columns: baseline_summary.loc[group, ('婚姻状况', f'{status}%')] = marriage_percentages.loc[group, status] if group in marriage_percentages.index else 0 for history in med_history_percentages.columns: baseline_summary.loc[group, ('既往用药', f'{history}%')] = med_history_percentages.loc[group, history] if group in med_history_percentages.index else 0 for level in depression_percentages.columns: baseline_summary.loc[group, ('抑郁程度', f'{level}%')] = depression_percentages.loc[group, level] if group in depression_percentages.index else 0 print("基线特征总结表：") print(baseline_summary) baseline_summary.to_excel('问题1结果/问题1_基线特征总结表.xlsx') # 统计检验结果总结表 age_anova = ols('年龄_岁 ~ C(组别)', data=combined_df).fit() age_anova_table = sm.stats.anova_lm(age_anova, typ=2) marriage_chi2, marriage_p, marriage_dof, _ = chi2_contingency(marriage_counts) med_history_chi2, med_history_p, med_history_dof, _ = chi2_contingency(med_history_counts) depression_chi2, depression_p, depression_dof, _ = chi2_contingency(depression_counts) stat_tests = pd.DataFrame({ '基线特征': ['年龄', '婚姻状况', '既往用药', '抑郁程度'], '统计方法': ['ANOVA', '卡方检验', '卡方检验', '卡方检验'], '检验统计量': [age_anova_table['F'][0], marriage_chi2, med_history_chi2, depression_chi2], 'p值': [age_anova_table['PR(>F)'][0], marriage_p, med_history_p, depression_p], '效应量': [eta_squared, cramer_v, cramer_v, cramer_v], '效应量类型': ['Eta Squared', 'Cramer\'s V', 'Cramer\'s V', 'Cramer\'s V'], '组间差异显著性': ['显著' if p < 0.05 else '不显著' for p in [age_anova_table['PR(>F)'][0], marriage_p, med_history_p, depression_p]] }) stat_tests['检验统计量'] = stat_tests['检验统计量'].map(lambda x: f"{x:.4f}") stat_tests['p值'] = stat_tests['p值'].map(lambda x: f"{x:.4f}") stat_tests['效应量'] = stat_tests['效应量'].map(lambda x: f"{x:.4f}") print("统计检验结果总结：") print(stat_tests) stat_tests.to_excel('问题1结果/问题1_统计检验结果总结.xlsx', index=False) 将他转为matlab代码，。同时保证所有功能不变

#!/usr/bin/env python # coding: utf-8 import pandas as pd import numpy as np import matplotlib as mpl import matplotlib.pyplot as plt import seaborn as sns # 导入seaborn库 from sklearn.preprocessing import MinMaxScaler, OneHotEncoder from sklearn.metrics import silhouette_score # 导入轮廓系数指标 from sklearn.cluster import KMeans # KMeans模块 ## 设置属性防止中文乱码 mpl.rcParams['font.sans-serif'] = [u'SimHei'] mpl.rcParams['axes.unicode_minus'] = False # 加载数据 raw_data = pd.read_csv(r'./ad_performance.csv') # 数据审查 print("数据框的前几行:") print(raw_data.head(2)) print("\n数据类型分布:") print(raw_data.info()) print("\n原始数据基本描述性信息:") print(raw_data.describe().round(2).T) # 查看缺失值情况 na_cols = raw_data.isnull().any(axis=0) print("\n每一列是否具有缺失值:") print(na_cols) print("\n具有缺失值的行总记录数:") print(raw_data.isnull().sum().sort_values(ascending=False)) # 变量之间的相关性分析 print("\n原始数据相关性信息:") corr = raw_data.corr(numeric_only=True).round(2) print(corr.T) # 热力图 plt.figure(figsize=(10, 8)) sns.heatmap(corr, cmap='coolwarm', annot=True, fmt=".2f") plt.title('变量相关性热力图') plt.show() # 删除平均停留时间列 raw_data2 = raw_data.drop(['平均停留时间'], axis=1) # 类别变量取值 cols = ["素材类型", "广告类型", "合作方式", "广告尺寸", "广告卖点"] for x in cols: data = raw_data2[x].unique() print(f"\n变量【{x}】的取值有：\n{data}") # 字符串分类独热编码处理 ohe_matrix1 = pd.get_dummies(raw_data2[cols]) print("\n独热编码后的矩阵头5行:") print(ohe_matrix1.head(5)) # 数据标准化 sacle_matrix = raw_data2.iloc[:, 1:7] # 获得要转换的矩阵 model_scaler = MinMaxScaler() # 建立MinMaxScaler模型对象 data_scaled = model_scaler.fit_transform(sacle_matrix) # MinMaxScaler标准化处理 print("\n标准化后的数据样本（前5行）:") print(pd.DataFrame(data_scaled.round(2), columns=sacle_matrix.columns).head(5)) # 合并所有维度 X = np.hstack((data_scaled, ohe_matrix1.values)) # 建立模型 score_list = list() # 用来存储每个K下模型的平局轮廓系数 silhouette_int = -1 # 初始化的平均轮廓系数阀值 best_kmeans = None cluster_labels_k = None for n_clusters in range(2, 8): # 遍历从2到7几个有限组 model_kmeans = KMeans(n_clusters=n_clusters, random_state=42) # 建立聚类模型对象 labels_tmp = model_kmeans.fit_predict(X) # 训练聚类模型 silhouette_tmp = silhouette_score(X, labels_tmp) # 得到每个K下的平均轮廓系数 if silhouette_tmp > silhouette_int: # 如果平均轮廓系数更高 best_k = n_clusters # 保存最好的K silhouette_int = silhouette_tmp # 保存平均轮廓得分 best_kmeans = model_kmeans # 保存模型实例对象 cluster_labels_k = labels_tmp # 保存聚类标签 score_list.append([n_clusters, silhouette_tmp]) # 将每次K及其得分追加到列表 print('\n{:*^60}'.format('K值对应的轮廓系数:')) print(np.array(score_list)) # 打印输出所有K下的详细得分 print(f'\n最优的K值是:{best_k} \n对应的轮廓系数是:{silhouette_int}') # 将原始数据与聚类标签整合 cluster_labels = pd.DataFrame(cluster_labels_k, columns=['clusters']) # 获得训练集下的标签信息 merge_data = pd.concat((raw_data2, cluster_labels), axis=1) # 将原始处理过的数据跟聚类标签整合 print("\n合并后的数据框前几行:") print(merge_data.head()) # 计算每个聚类类别下的样本量和样本占比 clustering_count = merge_data.groupby('clusters').size().reset_index(name='counts') # 计算每个聚类类别的样本量 clustering_ratio = clustering_count.copy() clustering_ratio['percentage'] = (clustering_ratio['counts'] / len(merge_data)).round(2) # 计算每个聚类类别的样本量占比 print("\n每个聚类类别的样本量:") print(clustering_count) print("\n#" * 30) print("\n每个聚类类别的样本量占比:") print(clustering_ratio) # 计算各个聚类类别内部最显著特征值 cluster_features = [] # 空列表，用于存储最终合并后的所有特征信息 for line in range(best_k): # 读取每个类索引 label_data = merge_data[merge_data['clusters'] == line] # 获得特定类的数据 part1_data = label_data.iloc[:, 1:7] # 获得数值型数据特征 part1_desc = part1_data.describe().round(3) # 得到数值型特征的描述性统计信息 merge_data1 = part1_desc.iloc[2, :] # 得到数值型特征的均值 part2_data = label_data.iloc[:, 7:-1] # 获得字符串型数据特征 part2_desc = part2_data.describe(include='all') # 获得字符串型数据特征的描述性统计信息 merge_data2 = part2_desc.iloc[2, :] # 获得字符串型数据特征的最频繁值 merge_line = pd.concat((merge_data1, merge_data2), axis=0) # 将数值型和字符串型典型特征沿行合并 cluster_features.append(merge_line) # 将每个类别下的数据特征追加到列表 # 输出完整的类别特征信息 cluster_pd = pd.DataFrame(cluster_features).T # 将列表转化为矩阵 print('\n{:*^60}'.format('每个类别主要的特征:')) all_cluster_set = pd.concat((clustering_count, clustering_ratio[['percentage']], cluster_pd), axis=1) # 将每个聚类类别的所有信息合并 print(all_cluster_set) # 曲线图 - 每个集群的日均UV变化趋势 plt.figure(figsize=(12, 6)) for i in range(best_k): cluster_data = merge_data[merge_data['clusters'] == i] plt.plot(cluster_data.index, cluster_data['日均UV'], label=f'Cluster {i+1}') plt.xlabel('Index') plt.ylabel('日均UV') plt.title('各聚类类别日均UV变化趋势') plt.legend() plt.grid(True) plt.show() # 直方图 - 各集群的订单转化率分布 plt.figure(figsize=(12, 6)) for i in range(best_k): cluster_data = merge_data[merge_data['clusters'] == i] plt.hist(cluster_data['订单转化率'], bins=10, alpha=0.5, label=f'Cluster {i+1}') plt.xlabel('订单转化率') plt.ylabel('Frequency') plt.title('各聚类类别订单转化率分布') plt.legend() plt.grid(True) plt.show() # 扇形图 - 各集群的样本占比 plt.figure(figsize=(8, 8)) plt.pie(clustering_ratio['percentage'], labels=[f'Cluster {i+1}' for i in range(best_k)], autopct='%1.1f%%', startangle=140) plt.title('各聚类类别样本占比') plt.axis('equal') # Equal aspect ratio ensures that pie is drawn as a circle. plt.show() # 雷达图 - 各集群的主要特征对比 num_sets = cluster_pd.iloc[:best_k, :].astype(np.float64) # 获取要展示的数据 labels = num_sets.index # 设置要展示的数据标签 colors = ['#FF9999','#66B3FF','#99FF99','#FFCC99', '#CCCCFF', '#FF99CC'] # 定义不同类别的颜色 angles = np.linspace(0, 2 * np.pi, len(labels), endpoint=False) # 计算各个区间的角度 angles = np.concatenate((angles, [angles[0]])) # 建立相同首尾字段以便于闭合 # 标准化数据 scaler = MinMaxScaler() num_sets_max_min = scaler.fit_transform(num_sets.T).T # 转置后标准化再转置回来 print("\n每个类别主要的数值特征（标准化后）:") print(num_sets_max_min) # 画图 fig = plt.figure(figsize=(10, 10)) # 建立画布 ax = fig.add_subplot(111, polar=True) # 增加子网格，注意polar参数 # 画雷达图 for i in range(len(num_sets_max_min)): # 循环每个类别 data_tmp = num_sets_max_min[i, :] # 获得对应类数据 data = np.concatenate((data_tmp, [data_tmp[0]])) # 建立相同首尾字段以便于闭合 ax.plot(angles, data, 'o-', c=colors[i], label=f"第{i+1}类渠道") # 画线 ax.fill(angles, data, alpha=0.25, color=colors[i]) # 设置图像显示格式 ax.set_thetagrids(angles[:-1] * 180 / np.pi, labels, fontproperties="SimHei") # 设置极坐标轴 ax.set_title("各聚类类别显著特征对比", pad=20, fontproperties="SimHei") # 设置标题放置 ax.set_rlim(0, 1.2) # 设置坐标轴尺度范围 ax.grid(True) plt.legend(loc="upper right", bbox_to_anchor=(1.2, 1.0), prop={'family': 'SimHei'}) plt.tight_layout(rect=[0, 0, 0.8, 1]) plt.show() # 数据结论 print("\n数据结论:") print("所有的渠道被分为{}各类别，每个类别的样本量分别为：{}\n对应占比分别为：{}".format( best_k, clustering_count['counts'].values.tolist(), clustering_ratio['percentage'].values.tolist()))