飞道的博客

Python 运筹优化4 遗传算法3 geatpy 实践

2021-05-14 13:04 203人阅读  评论(0)

说明

上篇介绍了基本的例子，这篇做一个实践。这里再推荐一篇文章，排版更好一些。

内容

1 业务问题

假设存在这样一个业务问题：

  • 1 在爬虫服务器上起k个同样的爬虫，爬取目标网站
  • 2 对于某个爬虫来说，每次爬取时需要使用一个ip发起请求
  • 3 某一时刻，爬虫服务器需要根据爬虫的数量向代理请求ip， ip有每分钟和每天最大数量两种限制

假设

  • 1 不同机器起的爬虫，处理速度一样
  • 2 目标主机的响应稳定
  • 3 同一个代理，ip的性能稳定

推论

  • 爬虫的速度𝑐𝑖取决于代理ip的质量
  • 稳定时间 𝑇𝑖 指每个ip的稳定可用时间
  • 使用上限 𝑁𝑖 指每个代理每天可提供的最大ip数量， 𝑛𝑖 是每分钟代理可提供的最大ip数量
  • 以平均每秒爬取企业数来衡量爬虫的速度 𝑐𝑖

目标

  • 优化代理池和爬虫服务器的数量，在一定成本下获得最高的爬取效率（爬取速度最快）

2 分析

假设两类服务器（高配和低配），机器1起 𝑘1 个爬虫，机器2起 𝑘2 个爬虫，每台服务器金额分别为 𝑚1 和 𝑚2 ，购置 𝑥1 和 𝑥2 台

假设存在两种代理，稳定时间分别为 𝑡1 和 𝑡2 ，每日使用上限分别为 𝑁1 和 𝑁2
单位时间可维持的最大ip数为 𝑛1 和 𝑛2 ，金额分别为 𝑀1 和 𝑀2 ，爬取速度为 𝑐1 和 𝑐2

分配给两种代理的爬虫数分别为 𝑦1 和 𝑦2
总预算为 𝑀
基于以上描述,可得：

优化变量：购买两类服务器（高配、低配）的数量 𝑥1,𝑥2 ，爬虫分配给两家代理的数量 𝑦1,𝑦2
参数：服务器起的爬虫数 𝑘1,𝑘2 ，服务器金额 𝑚1,𝑚2 ，两种代理稳定时间 𝑡1,𝑡2 ，每日使用上限 𝑁1,𝑁2 ，单位时间可维持的最大ip数 𝑛1,𝑛2 ，代理金额 𝑀1,𝑀2 ，爬虫爬取速度 𝑐1,𝑐2 ，总预算 𝑀

3 整数规划

之前用整数规划的方法完成过，看起来比较数学

from pulp import *
# 函数定义
# 定义函数，输入各项参数，返回变量 x1, x2, y1, y2取值和最优值

#参数：
#服务器起的爬虫数 k1,k2，服务器金额 m1, m2
#两种代理稳定时间 t1, t2，每日使用上限 N1, N2，单位时间可维持的最大ip数 n1, n2，代理金额 M1, M2
#爬虫爬取速度 c1, c2，总预算 M
def linprog(k1, k2, m1, m2, t1, t2, N1, N2, n1, n2, M1, M2, c1, c2, M):
   
    # Create a new model
    f = LpProblem(name="Problem1", sense=LpMaximize)

    # Create variables
    y1 = LpVariable(cat=LpInteger, name="y1", lowBound=0)
    y2 = LpVariable(cat=LpInteger, name="y2", lowBound=0)
    x1 = LpVariable(cat=LpInteger, name="x1", lowBound=0)
    x2 = LpVariable(cat=LpInteger, name="x2", lowBound=0)

    # Set objective
    f += c1 * y1 + c2 * y2, 'Obj'

    # Add constraint
    f += y1 + y2 - k1*x1 - k2*x2 == 0, 'c0'
    f += m1 * x1 + m2 * x2 <= M - M1 - M2 , 'c1'
    f += 24*60/t1 * y1 <= N1, 'c2'
    f += 24*60/t2 * y2 <= N2, 'c3'
    f += y1 <= n1, 'c4'
    f += y2 <= n2, 'c5'

    # Calculate
    status = f.solve()
    
    res = {
   }
    
    for v in f.variables():
        res[v.name] = v.varValue
    
    res['objective'] = value(f.objective)
    
    return res

接下来给参数赋值

#服务器起的爬虫数 k1,k2
k1 = 4
k2 = 3

#服务器金额 m1, m2
m1 = 2000
m2 = 600

#两种代理稳定时间 t1, t2
t1 = 2
t2 = 6

#每日使用上限 N1, N2
N1 = 150000
N2 = 580000

#单位时间可维持的最大ip数 n1, n2
n1 = 10
n2 = 10

#代理金额 M1, M2
M1 = 1000
M2 = 1360

#爬虫爬取速度 c1, c2
c1 = 1.2033
c2 = 1.9567

#总预算 M
M = 8000

执行函数并产生结果

# 测下时间 %%timeit
res = linprog(k1, k2, m1, m2, t1, t2, N1, N2, n1, n2, M1, M2, c1, c2, M)
---
16.1 ms ± 79.9 µs per loop (mean ± std. dev. of 7 runs, 100 loops each

print('x1 高配爬虫服务器（4个爬虫）:%s \nx2 低配爬虫服务器（3个爬虫） %s \ny1 分配给代理商1的爬虫 %s \ny2 分配给代理商2的爬虫 %s \n每秒爬取量：%.2f(%.2f/每天) \n' 
      %(res['x1'], res['x2'], res['y1'], res['y2'], res['objective'],res['objective']*86400))

x1 高配爬虫服务器（4个爬虫）:1.0 
x2 低配爬虫服务器（3个爬虫） 5.0 
y1 分配给代理商1的爬虫 9.0 
y2 分配给代理商2的爬虫 10.0 
每秒爬取量：30.40(2626274.88/每天)

速度还是相当快的，16ms。

4 遗传算法

此时的染色体x1,x2,x3,x4对应于整数规划里的x1,x2,y1,y2

import numpy as np
import geatpy as ea # 导入geatpy库
import time

#服务器起的爬虫数 k1,k2
k1 = 4
k2 = 3

#服务器金额 m1, m2
m1 = 2000
m2 = 600

#两种代理稳定时间 t1, t2
t1 = 2
t2 = 6

#每日使用上限 N1, N2
N1 = 150000
N2 = 580000

#单位时间可维持的最大ip数 n1, n2
n1 = 10
n2 = 10

#代理金额 M1, M2
M1 = 1000
M2 = 1360

#爬虫爬取速度 c1, c2
c1 = 1.2033
c2 = 1.9567

#总预算 M
M = 8000

def aimfunc(Phen, CV):
    x1 = Phen[:, [0]] 
    x2 = Phen[:, [1]]
    x3 = Phen[:, [2]]
    x4 = Phen[:, [3]]
#     x1, x2 , x3 , x4 = Phen.T
    f = c1* x3 + c2 * x4
    
    # 约束
    # 约束1：k1x1 + k2x2 - y1 - y2 = 0
    cv1 = np.abs(k1* x1 + k2* x2 - x3 - x4)
    # 约束2： x1m1 + x2m2 + M1 + M2 - M <=0
    cv2 = x1 * m1 + x2 * m2 + M1 + M2 - M
    # 约束3：x3 * 1day/ t1 - N1 
    cv3 = 1440 * x3/t1 - N1
    # 约束4：x4 * 1day/t2 - N2
    cv4 = 1440 * x4 / t2 - N2
    # 约束5：y1 -n1 <=0
    cv5 = x3 - n1
    # 约束6：y2 -n2 <=0
    cv6 = x4 - n2
    
    CV = np.hstack([cv1,cv2,cv3,cv4,cv5,cv6])
    return [f, CV]

# 变量设置
x1 = [0, 50] # 第一个决策变量范围
x2 = [0, 50] # 第二个决策变量范围
x3 = [0, 50]
x4 = [0, 50]


b1 = [1, 1] # 第一个决策变量边界，1表示包含范围的边界，0表示不包含
b2 = [1, 1] # 第二个决策变量边界，1表示包含范围的边界，0表示不包含
b3 = [1, 1]
b4 = [1, 1]


ranges=np.vstack([x1, x2,x3,x4]).T # 生成自变量的范围矩阵，使得第一行为所有决策变量的下界，第二行为上界
borders=np.vstack([b1, b2,b3,b4]).T # 生成自变量的边界矩阵
varTypes = np.array([1,1,1,1]) # 决策变量的类型，0表示连续，1表示离散

# 染色体编码设置
Encoding = 'BG' # 'BG'表示采用二进制/格雷编码
codes = [0, 0, 0,0] # 决策变量的编码方式，设置两个0表示两个决策变量均使用二进制编码
precisions =[4, 4,4,4] # 决策变量的编码精度，表示二进制编码串解码后能表示的决策变量的精度可达到小数点后6位
scales = [0, 0,0,0] # 0表示采用算术刻度，1表示采用对数刻度
FieldD = ea.crtfld(Encoding,varTypes,ranges,borders,precisions,codes,scales) # 调用函数创建译码矩阵

# 遗传算法参数设置
NIND      = 1000; # 种群个体数目
MAXGEN    = 2000; # 最大遗传代数
maxormins = [-1] # 列表元素为1则表示对应的目标函数是最小化，元素为-1则表示对应的目标函数是最大化
selectStyle = 'rws' # 采用轮盘赌选择
recStyle  = 'xovdp' # 采用两点交叉
mutStyle  = 'mutbin' # 采用二进制染色体的变异算子
pc        = 0.7 # 交叉概率
pm        = 1 # 整条染色体的变异概率（每一位的变异概率=pm/染色体长度）
Lind = int(np.sum(FieldD[0, :])) # 计算染色体长度
obj_trace = np.zeros((MAXGEN, 2)) # 定义目标函数值记录器
var_trace = np.zeros((MAXGEN, Lind)) # 染色体记录器，记录历代最优个体的染色体

start_time = time.time() # 开始计时
Chrom = ea.crtpc(Encoding, NIND, FieldD) # 生成种群染色体矩阵
variable = ea.bs2ri(Chrom, FieldD) # 对初始种群进行解码
CV = np.zeros((NIND, 1)) # 初始化一个CV矩阵（此时因为未确定个体是否满足约束条件，因此初始化元素为0，暂认为所有个体是可行解个体）
ObjV, CV = aimfunc(variable, CV) # 计算初始种群个体的目标函数值
FitnV = ea.ranking(maxormins * ObjV, CV) # 根据目标函数大小分配适应度值
best_ind = np.argmax(FitnV) # 计算当代最优个体的序号
# 开始进化
for gen in range(MAXGEN):
    SelCh = Chrom[ea.selecting(selectStyle,FitnV,NIND-1),:] # 选择
    SelCh = ea.recombin(recStyle, SelCh, pc) # 重组
    SelCh = ea.mutate(mutStyle, Encoding, SelCh, pm) # 变异
    # 把父代精英个体与子代的染色体进行合并，得到新一代种群
    Chrom = np.vstack([Chrom[best_ind, :].astype(int), SelCh])
    Phen = ea.bs2ri(Chrom, FieldD) # 对种群进行解码(二进制转十进制)
    ObjV, CV = aimfunc(Phen, CV) # 求种群个体的目标函数值
    FitnV = ea.ranking(maxormins * ObjV, CV) # 根据目标函数大小分配适应度值
    # 记录
    best_ind = np.argmax(FitnV) # 计算当代最优个体的序号
    obj_trace[gen,0]=np.sum(ObjV)/ObjV.shape[0] #记录当代种群的目标函数均值
    obj_trace[gen,1]=ObjV[best_ind] #记录当代种群最优个体目标函数值
    var_trace[gen,:]=Chrom[best_ind,:] #记录当代种群最优个体的染色体
# 进化完成
end_time = time.time() # 结束计时
# ea.trcplot(obj_trace, [['种群个体平均目标函数值', '种群最优个体目标函数值']]) # 绘制图像
ea.trcplot(obj_trace, [['Avg Target', 'Best Target']]) # 绘制图像


# 结果
best_gen = np.argmax(obj_trace[:, [1]])
print('最优解的目标函数值：', obj_trace[best_gen, 1])
variable = ea.bs2ri(var_trace[[best_gen], :], FieldD) # 解码得到表现型（即对应的决策变量值）
print('最优解的决策变量值为：')
for i in range(variable.shape[1]):
    print('x'+str(i)+'=',variable[0, i])
print('用时：', end_time - start_time, '秒')
---
最优解的目标函数值： 30.396700000000003
最优解的决策变量值为：
x0= 1
x1= 5
x2= 9
x3= 10
用时： 2.0616087913513184 秒

5 总结

  • 1 遗传算法和整数规划产生了相同的结果
  • 2 模拟的种群数量或者迭代次数很少会出现很大的偏差，甚至不能满足等式的约束
  • 3 限定x1,x2,x3,x4的返回，可以适当减少进化的次数（也能达到最优结果）
  • 4 实际使用时可以用对象封装，不会那么多代码
  • 5 为什么要用遗传算法？（比整数规划方法慢很多）
    • 1 采用Geatpy提供的进化算法框架可以既能最大程度地描述清楚所要求解的问题，而且与进化算法是高度脱耦的，即上面在编写问题类的时候完全不需要管后面采用什么算法、采用什么样编码的种群，只需把问题描述清楚即可。
    • 2 而且，遗传算法有个好处是：目标函数可以写得相当复杂，可以解决各种复杂的问题，比如神经网络。以BP神经网络为例，可以把神经网络的参数作为决策变量，神经网络的训练误差作为目标函数值，只需把上面的例子修改一下就行了

其实根本的原因在于，整数规划属于解析方法，而遗传算法属于模拟算法。模拟算法通用性更好，只是多耗费一些算力罢了。类似于深度学习和传统机器学习的关系。


转载：https://blog.csdn.net/yukai08008/article/details/116770162
查看评论
* 以上用户言论只代表其个人观点，不代表本网站的观点或立场