【Python】读取r语言数据+NMF算法（完整代码+详细注释）

最新推荐文章于 2023-10-25 11:28:01 发布

发现你走远了

最新推荐文章于 2023-10-25 11:28:01 发布

阅读量1.1k

点赞数 1

分类专栏： python 文章标签： python r语言算法 nmf

本文链接：https://blog.csdn.net/u011027547/article/details/130280801

版权

python 专栏收录该内容

201 篇文章 155 订阅

订阅专栏

依赖库

import numpy as np
import torch
import random
import matplotlib.pyplot as plt

import pyreadr
import pandas as pd

代码功能

导入r语言radata数据
通过修改切片选择数据区域，得到待分解的X矩阵
调用nmf函数时传入r, maxiter, minError

# 封装好的NMF算法
# X 传入的要被分解的矩阵
# r U V 的维度
# maxiter 最大迭代次数
# minError 迭代损失
def nmf(X, r, maxiter, minError):

数据持久化保存
绘制损失图

完整代码

# @Time    : 2023/4/19 8:19
# @Author  : 南黎
# @FileName: 基于pyreadr读入r语言数据转为csv数据.py

import numpy as np
import torch
import random
import matplotlib.pyplot as plt


# 封装好的NMF算法
# X 传入的要被分解的矩阵
# r U V 的维度
# maxiter 最大迭代次数
# minError 迭代损失
def nmf(X, r, maxiter, minError):
    # X=U*V'
    row, col = X.shape
    U = np.around(np.array(np.random.rand(row, r)), 5)#随机生成的U矩阵 服从“0~1”均匀分布的随机样本值。随机样本取值范围是[0,1)
    V = np.around(np.array(np.random.rand(col, r)), 5)#随机生成的V矩阵 服从“0~1”均匀分布的随机样本值。随机样本取值范围是[0,1)
    obj = []
    # 通过迭代
    for iter in range(maxiter):
        print('-----------------------------')
        print('开始第', iter, '次迭代')
        # update U
        XV = np.dot(X, V)# np.dot(a ,b), 其中a和b都是二维矩阵，此时dot就是进行的矩阵乘法运算
        UVV = np.dot(U, np.dot(V.T, V))
        U = (U * (XV / np.maximum(UVV, 1e-10)))
        # update V
        XU = np.dot(X.T, U)
        VUU = np.dot(V, np.dot(V.T, V))
        V = (V * (XU / np.maximum(VUU, 1e-10)))
        d = np.diag(1 / np.maximum(np.sqrt(np.sum(V * V, 0)), 1e-10))
        V = np.dot(V, d)

        temp = X - np.dot(U, np.transpose(V))#计算损失
        error = np.sum(temp * temp)#损失和
        print('error:', error)
        print('第', iter, '次迭代结束')
        obj.append(error)
        if error < minError:
            break

    return U, V, obj


import pyreadr
import pandas as pd
in_file_path="TCGA_BRCA_expr_raw.RData"#输入文件名
datas = pyreadr.read_r(in_file_path)
#print(type(datas))# <class 'collections.OrderedDict'>
# print(datas.keys()) #输出数据名 odict_keys(['BRCA.expr']) 本文的数据只有 BRCA.expr 有些文件可能是一个列表 有多个数据项
df = datas["BRCA.expr"] # 根据数据名得到数据，数据类型是 <class 'pandas.core.frame.DataFrame'>
# print(type(df))# <class 'pandas.core.frame.DataFrame'>
# print(df) #打印数据展示

# #展示列名
# col_names = df.columns
# print("展示列名\n",col_names)
#
# #展示行名
# index_names = df.index
# print("展示行名\n",index_names)

# # 取前面X行数据
# test_df=df.head(10)#取前面10行数据
# print("取前面10行数据\n",test_df)#取前面10行数据

# # df.iloc方法按照切片取指定X行数据
# # 注意：pandas的1.0.0版本后，已经对ix函数进行了升级和重构。   老版本这里不是df.iloc而是df.ix
# # 现在都是新版本了，如果你看到 df.ix这种操作也不要惊讶


#设置要分解矩阵的大小
test_df=df.iloc[0:20, 0:50]  # 使用 iloc 函数 ，切片类似matlab  取得第0-20行，第0-50列 （左闭右开）
print("df.iloc方法按照切片取指定X行数据\n",test_df)#取指定X行数据   [5 rows x 1222 columns]

test_df=test_df.values# [20 rows x 50 columns]
print(type(test_df)) # <class 'numpy.ndarray'>

X = np.array(np.abs(test_df))    #  确保非负
print("X shape\n", X.shape)# (20, 50)
# print('X:',X)#初始待分解的矩阵

U, V, obj = nmf(X, 2, 100, 0.01)#要分解的矩阵X； U, V的另一个维度 为2； 迭代次数100次  ；迭代损失0.01；
print("U\n",U)#     # X=U*V'
print("U shape\n", U.shape)#(20, 2)
out_file_path="U.csv"
print("开始转化"+out_file_path)
pd.DataFrame(U).to_csv(out_file_path)
print("转化完成"+out_file_path)#这个转化过程比较久，没有出现输出转化完成之前不要动文件，可能因为文件的只读导致文件损坏！！

V=V.T
print("V\n",V)#     # X=U*V'   V矩阵好像是要转置的（如果我的数学没有错的话，这里大家审查一下，哈哈）
print("V shape\n", V.shape)# (50, 2)
out_file_path="V.csv"
print("开始转化"+out_file_path)
pd.DataFrame(V).to_csv(out_file_path)
print("转化完成"+out_file_path)#这个转化过程比较久，没有出现输出转化完成之前不要动文件，可能因为文件的只读导致文件损坏！！


x = range(len(obj))
# 绘图部分 随着迭代次数损失降低曲线
plt.plot(x, obj)
plt.show()