网络毁伤最大化问题

本文分享了如何使用Matlab和Python代码来实现通信网络的毁伤最大化问题,通过目标函数计算网络效率下降,对比了CELF、GABIN和度/接近度/介数排序算法。Matlab初试代码和Python高效版本的详细步骤在文中提供,适合研究者和开发者参考。
摘要由CSDN通过智能技术生成

通信网络毁伤最大化问题


论文已中,把刚开始写的代码放一下(最终版就不放了),版权所有,禁止抄袭,转载请注明出处。

目标函数

网络效率为 η = 1 N ( N − 1 ) ∑ v i ≠ v j ∈ V 1 d i j \eta =\frac{1}{{N(N - 1)}}\sum\limits_{v_i \ne v_j \in V}\frac{1}{{d_{ij}}} η=N(N1)1vi=vjVdij1

则衡量网络毁伤的目标函数即网络下降效率为
R ( s ) = ( 1 − η ( G ′ ) / η ( G ) ) × 100 % R(s) = (1-\eta (G')/\eta(G))\times100\% R(s)=(1η(G)/η(G))×100%

下面对比了基于拓扑势的CELF算法(TPCELF),基于重要节点的贪心算法(GABIN)和基于度值(Degree)、接近度(Closeness)和介数(Betweenness)排序的算法,Matlab程序是初步尝试,犹豫在计算介数的时候总存在一些时间消耗过大的问题,按照常用的计算介数方法编程消耗太大,对大型网络等候时间过长,Lev的工具箱暂时不好用,总是报错。这里改用了matlab_bgl工具箱,具体的地址见这里,使用的时候需要把private文件夹重命名一下在添加到路径。Matlab的实验比较简单,基于拓扑势的方法也没有加CELF。

后来发现python的好几个包在这块计算上做的比较好,改用基于python编程。

代码总体是完整的,只是没有放生成仿真数据的的内容,在另一个文档里,想比读者肯定能够自己实现的,非常easy,网上也有很多现成的数据集,改一下路径就能跑~

Matlab代码

close all
clear 
tic

%% 生成有无标度网络
%方式一 自己生成
% BA=Powerlaw_with_Expected_Exponent(300,6,2,1);%这个函数生成的幂律分布的网络会有孤立点存在

%方式二 用snap网络生成的数据集,如果是批量的,加一个循环就可以
m = 1;
filename = strcat('D:\test_300_1.5_',num2str(m),'.txt');
d=textread(filename,'','headerlines',3);%跳过第3行开始读数据
s=[];t=[];
for i=1:length(d)
    s=[s,d(i,1)+1];
    t=[t,d(i,2)+1];
end
g = graph(s,t);
BA = adjacency(g); %下文介数是要用稀疏矩阵的

A = BA;
A_1 = BA;                 %用于比较,保存初始矩阵,liu的方法

G =graph(BA);
% BA = adjacency(g);    %得到邻接矩阵
D = distances(G);      %得到距离矩阵
plot(G)                %展示图片是怎样的

%% 初始化相应的值
N=height(G.Nodes);    %节点数
DD = D+1e6*eye(N);    %在对角线上为了排除自己和自己的距离,令其为一个很大的数
Initial_efficiency = sum(sum(1./DD));    %初始的网络效率
S=[];   S_1=[];       %两个种子集
sigma = 0.9428;       %影响力因子,此处设置为两跳的最小值 
alpha = 8;            %Liu论文的节点初始集合放大系数
%% 想要打击K个节点

K=5;
for k=K: -1 :1
   
    %计算节点的拓扑势
    for i=1:N
        phi(i)=sum(exp(-(D(i,:)/(sigma)).^2));
    end
    [BB,I] = sort( phi,'descend');     %对拓扑势进行排序
    

    [B,E] = betweenness_centrality(A_1);        %计算介数
    Du = Degree_undirect(A_1);          %计算度数
    De = Du+max(B)/max(Du).*B;          %计算节点重要度
    [BB_1,I_1] = sort( De ,'descend');  %对重要度进行排序
    
    %这里取用最小alpha=10,相对固定住
    for i = 1:max(k,10)
        a=A;
        a(I(i),:)=0;
        a(:,I(i))=0;
        g = graph(a);
        D = distances(g);
        DD = D+1e6*eye(N);
        descend_efficiency(i) = 1-sum(sum(1./DD))/Initial_efficiency;
    end
    
    for i=1:alpha*k
        a_1=A_1;
        a_1(I(i),:)=0;
        a_1(:,I(i))=0;
        g_1 = graph(a_1);
        D_1 = distances(g_1);
        DD_1 = D_1+1e6*eye(N);
        descend_efficiency_1(i) = 1-sum(sum(1./DD_1))/Initial_efficiency;
    end
    
    ddd(K+1-k)=max(descend_efficiency);
    index=I(find(descend_efficiency == max(descend_efficiency),1));%找到最大衰减性能对应的节点的位置
    S = [S,index];  %把节点加入种子集合
    A(index,:)=0;   %删除被打击的节点
    A(:,index)=0;
    
    
    
    ddd_1(K+1-k)=max(descend_efficiency_1);
    index_1=I_1(find(descend_efficiency_1 == max(descend_efficiency_1),1));%找到最大衰减性能对应的节点的位置
    S_1 = [S_1,index_1];  %把节点加入种子集合
    A_1(index_1,:)=0;   %删除被打击的节点
    A_1(:,index_1)=0;

end
toc
figure
plot(ddd,'-*');hold on ;
plot(ddd_1)
legend('我','Liu')

Python代码

import networkx as nx
import matplotlib.pyplot as plt
import time
import numpy as np
import copy
import pandas as pd
import os

'''
=====计算增益的等价说明======
在计算每个节点删除后的增益时,采用了选择最大下降效率的方法,其实这和最大增益是等价的
因为最大下降效率 - 上一状态网络的下降效率 = 最大增益,被减数在每轮中是一个定值
'''

def txt_read(data):
    '''
    读取txt文档,生成网络
    :param data:
    :return:
    '''
    return nx.read_edgelist(data, create_using=nx.Graph())

def effectiveness(path):
    '''
    读取最短路径,计算网络的效率count
    :param path:
    :return:
    '''
    count = 0
    for x in path.keys():
        for y in path[x].keys():
            try:
                count += 1/path[x][y]
            except:
                pass
    return count

def CELF_algrithm(path, Use):
    '''
    CELF算法封装
    传入总体参数 g N
    Initial_efficiency
    :param path: 初始的最短路径
    :param Use: bool 是否用celf思想
    :return: S_1  decrease_efficiency_1  time1
    '''

    # ##### #####开始进行网络毁伤实验  Part1  ##### #####
    start1 = time.time()
    decrease_efficiency_1 = []
    S_1 = []

    k = 0
    while k < K:
        phi = []
        for i in range(N):  # 这样的用法是为了让节点按照顺序排列
            temp = 0
            for j in g.nodes():
                try:
                    temp += np.exp(-(path[j][str(i)] / (sigma)) ** 2)
                except KeyError:
                    # 有的节点之间没有最短路径,会报错,采用异常处理模式,这样就不影响结果
                    pass
            phi.append(temp)
        idx = sorted(enumerate(phi), key=lambda x: x[1], reverse=True)  # 降序排列,得到id的元组

        efficiency = []
        for p in range(max(K - k, flex_alpha)):
            a = copy.deepcopy(G_1)
            a.remove_node(str(idx[p][0]))  # 删除节点以及对应的边
            a.add_node(str(idx[p][0]))  # 添加被删除的节点
            # 得到网络的最短路径,保存为两级的字典格式
            path = dict(nx.shortest_path_length(a))
            # 计算网络下降效率
            efficiency.append(1 - effectiveness(path) / Initial_efficiency)

        id_top3 = pd.Series(efficiency).sort_values(ascending=False).index[:3]  # 获取效率排名前三的效率值的id

        decrease_efficiency_1.append(efficiency[id_top3[0]])  # 选择效率下降最大的节点,也可用max(efficiency)
        index = str(idx[efficiency.index(efficiency[id_top3[0]])][0])  # 下标是字符串类型的,效率下降对应的p值,再对应返回过去idx
        S_1.append(index)
        G_1.remove_node(index)  # 删除节点以及对应的边
        G_1.add_node(index)  # 添加被删除的节点

        # plt.figure()
        # nx.draw_networkx(G_1)
        if Use == True:

            # 采用 CELF 算法思想:B的新一轮的增益要大于C在本轮的增益
            # 判断一下是否已经到达零界位置,否则在k刚好满足要求的时候会多运行一步
            if k + 1 == K:
                break

            # 计算C在本轮的增益
            if k == 0:  # 如果是打击第一个节点
                delta_C = efficiency[id_top3[2]]
                # print('delta_C=', delta_C)
            else:
                delta_C = efficiency[id_top3[2]] - decrease_efficiency_1[k - 1]
                # print('delta_C=',delta_C)

            index_B_new = id_top3[1]
            a = copy.deepcopy(G_1)
            a.remove_node(str(idx[index_B_new][0]))  # 删除节点以及对应的边
            a.add_node(str(idx[index_B_new][0]))  # 添加被删除的节点
            # 得到网络的最短路径,保存为两级的字典格式
            path = dict(nx.shortest_path_length(a))
            # 计算网络下降效率增益

            delta_B_new = (1 - effectiveness(path) / Initial_efficiency) - decrease_efficiency_1[k]

            if delta_C < delta_B_new:
                # 如果新的一轮B的增益大于本轮C的增益,直接用B,然后删除节点
                decrease_efficiency_1.append(1 - effectiveness(path) / Initial_efficiency)
                S_1.append(str(idx[index_B_new][0]))
                G_1.remove_node(str(idx[index_B_new][0]))  # 删除节点以及对应的边
                G_1.add_node(str(idx[index_B_new][0]))  # 添加被删除的节点
                k += 2
                # plt.figure()
                # nx.draw_networkx(G_1)
            else:  #如果B的新增益比原来的C要小
                k += 1

        else: #如果不同celf思想
            k += 1

    end1 = time.time()
    time1 = end1 - start1

    return S_1, decrease_efficiency_1, time1

def GABIN_algrithm():
    '''

    :return:
    '''

    # ##### #####开始进行网络毁伤实验  Part2 ##### #####
    start2 = time.time()
    decrease_efficiency_2 = []
    S_2 = []
    k = 0
    while k < K:
        
        # 得到网络的介数,存储是字典
        B = nx.betweenness_centrality(G_2, normalized=False)
        # 得到的度值,存储是字典
        Du = nx.degree(G_2)
        # 计算最大介数和最大度之间的比值,重要度
        De = max(B.values()) / max(Du, key=lambda x: x[-1])[-1]
        # 计算节点重要度
        Importance_D_B = {}
        for key in B.keys():
            Importance_D_B[key] = Du[key] + De * B[key]
        # [Importance_D_B[key]] = [Du[key] + De * B.get(key) for key in B.keys()]  # 有问题,这样就把原来的序号给搞没了
        idx = sorted(Importance_D_B.items(), key=lambda x: x[1], reverse=True)  # 降序排列

        efficiency = []
        for p in range(alpha * (K - k)):
            a = copy.deepcopy(G_2)
            a.remove_node(str(idx[p][0]))  # 删除节点以及对应的边
            a.add_node(str(idx[p][0]))  # 添加被删除的节点
            # 得到网络的最短路径,保存为两级的字典格式
            path = dict(nx.shortest_path_length(a))
            # 计算网络下降效率
            efficiency.append(1 - effectiveness(path) / Initial_efficiency)

        decrease_efficiency_2.append(max(efficiency))  # 选择效率下降最大的节点
        index = str(idx[efficiency.index(max(efficiency))][0])  # 下标是字符串类型的,效率下降对应的p值,再对应返回过去idx
        S_2.append(index)
        G_2.remove_node(index)  # 删除节点以及对应的边
        G_2.add_node(index)  # 添加被删除的节点

        k += 1

    end2 = time.time()
    time2 = end2 - start2
    return S_2, decrease_efficiency_2, time2

def Greedy_algrithm():
    '''
    贪婪算法
    :param g:
    :return:
    '''

    # ##### #####开始进行网络毁伤实验  Part3 ##### #####
    start3 = time.time()
    decrease_efficiency_3 = []
    S_3 = []

    k = 0
    while k < K:
        efficiency = []
        for p in G_3.nodes():
            a = copy.deepcopy(G_3)
            a.remove_node(p)  # 删除节点以及对应的边
            a.add_node(p)  # 添加被删除的节点
            # 得到网络的最短路径,保存为两级的字典格式
            path = dict(nx.shortest_path_length(a))
            # 计算网络下降效率,保存成节点编号和效率值的tuple
            efficiency.append((p, 1 - effectiveness(path) / Initial_efficiency))

        decrease_efficiency_3.append(max(efficiency,key= lambda x:x[-1])[-1])  # 选择效率下降最大的节点
        index = max(efficiency,key= lambda x:x[-1])[0] # 下标是字符串类型的

        S_3.append(index)
        G_3.remove_node(index)  # 删除节点以及对应的边
        G_3.add_node(index)  # 添加被删除的节点

        k += 1

    end3 = time.time()
    time3 = end3 - start3
    return S_3, decrease_efficiency_3, time3

def Du_algrithm():
    '''

    :return:
    '''
    start4 = time.time()
    G_4 = copy.deepcopy(g)
    Du = nx.degree(G_4)
    temp = sorted(Du, key=lambda x: x[1], reverse=True)[0:K]  # 元组降序排列
    S_4 = [temp[key][0] for key in range(K)]   # 得到种子集

    decrease_efficiency_4 = []
    for p in S_4:
        G_4.remove_node(p)
        G_4.add_node(p)
        path = dict(nx.shortest_path_length(G_4))
        decrease_efficiency_4.append(1 - effectiveness(path) / Initial_efficiency)

    end4 = time.time()

def Closness_algrithm():
    '''

    :return:
    '''
    start5 = time.time()
    G_5 = copy.deepcopy(g)
    Closeness = nx.closeness_centrality(G_5)
    temp = sorted(Closeness.items(), key=lambda x: x[1], reverse=True)[0:K]  # 字典降序排列
    S_5 = [temp[key][0] for key in range(K)]  # 得到种子集

    decrease_efficiency_5 = []
    for p in S_5:
        G_5.remove_node(p)
        G_5.add_node(p)
        path = dict(nx.shortest_path_length(G_5))
        decrease_efficiency_5.append(1 - effectiveness(path) / Initial_efficiency)

    end5 = time.time()
    return S_5, decrease_efficiency_5

def Betweenness_algrithm():
    '''

    :return:
    '''
    start6 =time.time()
    G_6 = copy.deepcopy(g)
    Betweenness = nx.betweenness_centrality(G_6)
    temp = sorted(Betweenness.items(), key=lambda x: x[1], reverse=True)[0:K]  # 字典降序排列
    S_6 = [temp[key][0] for key in range(K)]  # 得到种子集
    decrease_efficiency_6 = []
    for p in S_6:
        G_6.remove_node(p)
        G_6.add_node(p)
        path = dict(nx.shortest_path_length(G_6))
        decrease_efficiency_6.append(1 - effectiveness(path) / Initial_efficiency)

    end6 = time.time()
    return S_6, decrease_efficiency_6


if __name__ == '__main__':

    '''
    实验1获取数据 设定gamma=2 K=5 alpha=5  flex_alpha=10
    ####################################################
    '''
    Num_of_Nodes = [50, 60, 70, 80, 90, 100, 200, 300, 400, 500, 600, 700]
    Gamma = 2
    fsave = open('result1.txt', mode='w+')
    for num in Num_of_Nodes:
        filename = 'D:\\PY\\TEST1\\test_' + str(num) + '_2.txt'
        g = txt_read(filename)
        '''
        ####################################################
        '''
        #     # 得到网络的最短路径,保存为两级的字典格式
        path = dict(nx.shortest_path_length(g))
        #     # 设置备份网络并展示原始网络,不能简单传递变量赋值,这是链接到内存地址的
        G_1 = copy.deepcopy(g)
        G_2 = copy.deepcopy(g)
        G_3 = copy.deepcopy(g)
    
        '''
        ##### ##### ##### ##### ##### ##### ##### #####
        初始化参数设置
        '''
        sigma = 0.9428
        # 节点个数
        N = g.number_of_nodes()
        # 需要毁伤的节点个数
        K = 5
        # 种子候选点的放大系数
        alpha = 5
        # 拓扑势方法单一选择的最大的值
        flex_alpha = 10
        # 初始的网络效率
        Initial_efficiency = effectiveness(path)
        '''
        ##### ##### ##### ##### ##### ##### ##### #####
        '''
    
         # 输出部分的内容
         out1 = CELF_algrithm(path, True)
         out2 = GABIN_algrithm()
         out3 = Greedy_algrithm() 

         fsave.write(str(num)+',  '+str(out1[2])+',  '+str(out2[2])+',  '+str(out3[2])+'\n')
    fsave.close()


    start = time.time()
    '''
    实验2 获取数据 设定gamma=2.5 N=300 K=1-10  alpha=5  flex_alpha=10
    ####################################################
    '''
    # 设置K
    K = 10
    fsave = open('K_=_10.txt', 'w+')
    for i in range(1, 101):
        filename = 'D:\\PY\\TEST2\\test_300_' + str(2.5) + '_' + str(i) + '.txt'
        g = txt_read(filename)

        # 得到网络的最短路径,保存为两级的字典格式
        path = dict(nx.shortest_path_length(g))

        # 设置备份网络并展示原始网络,不能简单传递变量赋值,这是链接到内存地址的
        G_1 = copy.deepcopy(g)
        G_2 = copy.deepcopy(g)

        '''
        ##### ##### ##### ##### ##### ##### ##### #####
        初始化参数设置
        '''
        sigma = 0.9428
        # 节点个数
        N = g.number_of_nodes()
        # 种子候选点的放大系数
        alpha = 5
        # 拓扑势方法单一选择的最大的值
        flex_alpha = 10
        # 初始的网络效率
        Initial_efficiency = effectiveness(path)
        '''
        ##### ##### ##### ##### ##### ##### ##### #####
        '''
        out1 = CELF_algrithm(path, True)
        out2 = GABIN_algrithm()
        out3 = Du_algrithm()
        out4 = Closness_algrithm()
        out5 = Betweenness_algrithm()

        fsave.write(str(out1) + '\n')
        fsave.write(str(out2) + '\n')
        fsave.write(str(out3) + '\n')
        fsave.write(str(out4) + '\n')
        fsave.write(str(out5) + '\n')
        fsave.write('\n')

    fsave.close()

    end = time.time()

    print(end-start)

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

Haleine

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值