当前位置:   article > 正文

马尔可夫随机场的python实现_scipy随机场

scipy随机场

       前几天看图网络的时候,不知怎么从贝叶斯网络突然跳到了马尔科夫随机场,感觉还有点意思,不过还是没有完全理清其中的逻辑,网上讲的比较乱,参考博主on2way这篇博文《从贝叶斯理论到图像马尔科夫随机场》我梳理了一下思路,找时间总结一下,先挖一个坑吧,由于博主给出的是matlab实现,我本身不怎么使用matlab,于是就简单的用python复现了一下思路,由于没有使用kmeans算法初始化,结果差距比较明显,但是还算不错,其中应该是有一些问题,包括最后多分类的接口没有调好,还有四分类显示结果貌似不太对,由于后天要考英语今天看来是没法调完了,依旧是放着这个坑,有时间的时候代码改一下,以下是python3的代码,直接用jupyter notebook写的:

  1. import numpy as np
  2. import cv2
  3. import copy
  4. from scipy.stats import norm
  5. import matplotlib as plot
  6. img = cv2.imread('Lena.jpeg',0)
  7. cluster_num = 4
  8. maxiter = 50
  9. #随机初始化标签
  10. label = np.random.randint(0,cluster_num,img.shape)
  11. def cal_pp(label):
  12. result = np.zeros((cluster_num, 8, label.shape[0], label.shape[1]))
  13. statis_result = np.zeros((cluster_num, label.shape[0], label.shape[1]))
  14. for i in range(cluster_num):
  15. label_temp = copy.deepcopy(label)
  16. #先设置为4再设置为0,防止处理0时出错
  17. label_temp[label_temp!=i] = 4
  18. label_temp[label_temp==i] = 1
  19. label_temp[label_temp==4] = 0
  20. #根据label_temp可以计算出8个方向的标签
  21. #从左上角起顺时针方向,依次为0-8
  22. result[i][0][:,0] = 0
  23. result[i][0][0,:] = 0
  24. result[i][0][1:,1:] = label_temp[:-1,:-1]
  25. result[i][1][:,0] = 0
  26. result[i][1][1:,:] = label_temp[:-1,:]
  27. result[i][2][:,0] = 0
  28. result[i][2][0,:] = 0
  29. result[i][2][1:,1:] = label_temp[1:,:-1]
  30. result[i][3][:,0] = 0
  31. result[i][3][:,1:] = label_temp[:,:-1]
  32. result[i][4][:,-1] = 0
  33. result[i][4][:,:-1] = label_temp[:,1:]
  34. result[i][5][:,0] = 0
  35. result[i][5][-1,:] = 0
  36. result[i][5][1:,:-1] = label_temp[1:,:-1]
  37. result[i][6][-1,:] = 0
  38. result[i][6][:-1,:] = label_temp[1:,:]
  39. result[i][7][:,-1] = 0
  40. result[i][7][-1,:] = 0
  41. result[i][7][:-1,:-1] = label_temp[1:,1:]
  42. statis_result[i] = result[i][0]+result[i][1]+result[i][2]+result[i][3]+result[i][4]+result[i][5]+result[i][6]+result[i][7]
  43. statis_result = np.array(statis_result,dtype=np.float) / 8
  44. #防止出现0,原因不清楚
  45. statis_result[statis_result==0] = 0.001
  46. return statis_result
  47. def cal_lf(img, label):
  48. result = np.zeros((cluster_num, 8, label.shape[0], label.shape[1]))
  49. distribution_parameter = np.zeros((cluster_num,2),np.float)
  50. for i in range(cluster_num):
  51. img_temp = copy.deepcopy(img)
  52. img_temp[label!=i] = 0
  53. #根据原图灰度统计每一个标签的分布和方差
  54. img_list = img_temp.tolist()
  55. img_list = [item for sublist in img_list for item in sublist if item!= 0]
  56. #再根据原图的灰度数值算出似然函数值
  57. distribution_parameter[i][0] = np.mean(img_list)
  58. distribution_parameter[i][1] = np.std(img_list)
  59. #计算每一点属于不同类型的概率
  60. result_lf = np.zeros((cluster_num,label.shape[0],label.shape[1]),np.float)
  61. for i in range(cluster_num):
  62. result_lf[i] = norm(distribution_parameter[i][0],distribution_parameter[i][1]).pdf(img)
  63. return result_lf
  64. def update_label(img, init_label, maxiter):
  65. label = init_label
  66. for i in range(maxiter):
  67. prior_probability = cal_pp(label)
  68. likelihood = cal_lf(img,label)
  69. probability = np.zeros((cluster_num,label.shape[0],label.shape[1]))
  70. for j in range(cluster_num):
  71. probability[j] = prior_probability[j]*likelihood[j]
  72. #根据概率值更新label
  73. max_probability = probability.max(0)
  74. #label_one不需要,这一部分代码需要做一下调整,时间紧迫先演示一下
  75. label_two = probability[1]-max_probability
  76. label_two[label_two!=0] = 4
  77. label_two[label_two==0] = 1
  78. label_two[label_two==4] = 0
  79. label_three = probability[2]-max_probability
  80. label_three[label_three!=0] = 4
  81. label_three[label_three==0] = 2
  82. label_three[label_three==4] = 0
  83. label_four = probability[3]-max_probability
  84. label_four[label_four!=0] = 4
  85. label_four[label_four==0] = 3
  86. label_four[label_four==4] = 0
  87. label= label_two + label_three + label_four
  88. #每次绘制一张图片
  89. cv2.imshow("test",label)
  90. cv2.waitKey(500)
  91. cv2.destroyAllWindows()
  92. return label
  93. #main函数中运行,博主使用jupyter notebook就没有改了,请需要的自行改一下
  94. final = update_label(img,label,maxiter)

原图放这里供大家下载使用:

贴一下大致效果:
其中前两排是第1-8次迭代结果,最后一排是第12、16、20、24次迭代结果,效果还行

声明:本文内容由网友自发贡献,不代表【wpsshop博客】立场,版权归原作者所有,本站不承担相应法律责任。如您发现有侵权的内容,请联系我们。转载请注明出处:https://www.wpsshop.cn/w/Guff_9hys/article/detail/915997
推荐阅读
相关标签
  

闽ICP备14008679号