Moo*_*per 7 optimization r cost-management loss-function
我有属于不同类别的个体,他们位于不同的区域,这些人口预计会从population
低于该demand
值的值增长到该值。
population_and_demand_by_category_and_zone <- tibble::tribble(
~category, ~zone, ~population, ~demand,
"A", 1, 115, 138,
"A", 2, 121, 145,
"A", 3, 112, 134,
"A", 4, 76, 91,
"B", 1, 70, 99,
"B", 2, 59, 83,
"B", 3, 86, 121,
"B", 4, 139, 196,
"C", 1, 142, 160,
"C", 2, 72, 81,
"C", 3, 29, 33,
"C", 4, 58, 66,
"D", 1, 22, 47,
"D", 2, 23, 49,
"D", 3, 16, 34,
"D", 4, 45, 96
)
Run Code Online (Sandbox Code Playgroud)
区域具有给定的容量,当前人口低于此阈值,但某些区域的需求将超过容量。
demand_and_capacity_by_zone <- tibble::tribble(
~zone, ~demand, ~capacity, ~capacity_exceeded,
1, 444, 465, FALSE,
2, 358, 393, FALSE,
3, 322, 500, FALSE,
4, 449, 331, TRUE
)
Run Code Online (Sandbox Code Playgroud)
所以我们需要将这些人转移到一个新的区域(我们假设我们有足够的总容量)。我们需要移动的每个人都会产生成本,这取决于其类别和目的地区域。这些费用如下。
costs <- tibble::tribble(
~category, ~zone, ~cost,
"A", 1, 0.1,
"A", 2, 0.1,
"A", 3, 0.1,
"A", 4, 1.3,
"B", 1, 16.2,
"B", 2, 38.1,
"B", 3, 1.5,
"B", 4, 0.1,
"C", 1, 0.1,
"C", 2, 12.7,
"C", 3, 97.7,
"C", 4, 46.3,
"D", 1, 25.3,
"D", 2, 7.7,
"D", 3, 67.3,
"D", 4, 0.1
)
Run Code Online (Sandbox Code Playgroud)
我希望找到跨区域和类别的个体分布,以便使总成本最小化。所以基本上new_population
在population_and_demand_by_category_and_zone
上面描述的表中有一个新列。
如果有几个可能的解决方案,任何一个都可以,如果结果是一个非整数总体就好了。
真正的用例有大约 40 个类别和区域,所以更大但不是那么大。
这似乎是一个很常见的问题,所以我希望在 R 中有一种方便的方法来解决这个问题。
谢谢你的帮助!
这可以建模为小型 LP(线性规划)模型。我们引入非负变量move(c,z,z')
来指示从区域 z 移动到区域 z' 的类别 c 的人数。数学模型可以如下所示:
这可以使用任何 LP 求解器来实现。解决方案可能如下所示:
---- 83 VARIABLE move.L moves needed to meet capacity
zone1 zone2 zone3
catA.zone1 6
catA.zone4 29 62
catC.zone4 27
---- 83 VARIABLE alloc.L new allocation
zone1 zone2 zone3 zone4
catA 132 180 196
catB 99 83 121 196
catC 187 81 33 39
catD 47 49 34 96
---- 83 VARIABLE totcost.L = 12.400 total cost
Run Code Online (Sandbox Code Playgroud)
笔记:
allocation = demand + inflow - outflow
move(c,z,z)=0
确保我们不会从 z 移动到 z 本身。这个约束并不是真正需要的(它是由成本隐式强制执行的)。为了清楚起见,我添加了它。实际上,我通过将 的上限设置move(c,z,z)
为零(即没有显式约束)来实现这一点。对于非常大的模型,我会使用另一种可能性:甚至不生成变量move(c,z,z)
。这个模型很小,所以没有必要。如果您愿意,您可以完全省略它。population
模型中使用。我认为没有必要,除非我们看看下一个项目符号。我采用了 Erwin 的公式,对其进行了修改,以考虑分配应该大于每个区域和类别的人口(这意味着已经存在的个体不会移动),并使用 {lpSolve} 包实现它,这不会'不需要安装外部系统库。
\n埃尔文的解决方案可以通过以下方式获得move_new_only <- FALSE
。
library(tidyverse)\nlibrary(lpSolve)\n\nmove_new_only <- TRUE # means population in place can't be reallocated\ncategories <- unique(population_and_demand_by_category_and_zone$category)\nzones <- unique(population_and_demand_by_category_and_zone$zone)\nn_cat <- length(categories)\nn_zones <- length(zones)\n\n# empty coefficient arrays\nmove_coefs_template <- array(0, c(n_zones, n_zones, n_cat), \n dimnames = list(zones, zones, categories))\nalloc_coefs_template <- matrix(0, n_zones, n_cat, \n dimnames = list(zones, categories))\n\nbuild_zone_by_category_matrix <- function(data, col) {\n data %>% \n pivot_wider(\n id_cols = zone, names_from = category, values_from = {{col}}) %>% \n as.data.frame() %>% \n `row.names<-`(.$zone) %>% \n select(-zone) %>% \n as.matrix()\n}\n\ndemand_mat <- build_zone_by_category_matrix(\n population_and_demand_by_category_and_zone, demand)\n\ncost_mat <- build_zone_by_category_matrix(costs, cost)\n\npopulation_mat <- build_zone_by_category_matrix(\n population_and_demand_by_category_and_zone, population)\n
Run Code Online (Sandbox Code Playgroud)\n# stack the cost matrix vertically to build an array of all move coefficients\ncoefs_obj <- move_coefs_template\nfor(i in 1:n_zones) {\n coefs_obj[i,,] <- cost_mat\n}\n# flatten it for `lp`s `objective.in` argument, adding alloc coefs\nf.obj <- c(coefs_obj, alloc_coefs_template)\n
Run Code Online (Sandbox Code Playgroud)\ncoefs_con <- list() \nfor (z in zones) {\n coefs_con_zone <- list() \n for(categ in categories) {\n coefs_arrivals <- move_coefs_template\n coefs_arrivals[,z, categ] <- 1\n coefs_departures <- move_coefs_template\n coefs_departures[z,, categ] <- 1\n coefs_moves <- coefs_arrivals - coefs_departures\n coefs_alloc <- alloc_coefs_template\n coefs_alloc[z, categ] <- -1\n # flatten the array\n coefs_con_zone[[categ]] <- c(coefs_moves, coefs_alloc)\n }\n coefs_con[[z]] <- do.call(rbind, coefs_con_zone)\n}\n\n# stack the flattened arrays to build a matrix\nf.con1 <- do.call(rbind, coefs_con)\nf.dir1 <- rep("==", n_zones * n_cat)\nf.rhs1 <- -c(t(demand_mat)) # transposing so we start with all zone 1 and so on\n
Run Code Online (Sandbox Code Playgroud)\ncoefs_con <- list() \nfor (z in zones) {\n coefs_alloc <- alloc_coefs_template\n coefs_alloc[z, ] <- 1\n coefs_con[[z]] <- c(move_coefs_template, coefs_alloc)\n}\n\n# stack the flattened arrays to build a matrix\nf.con2 <- do.call(rbind, coefs_con)\nf.dir2 <- rep("<=", n_zones)\n\nf.rhs2 <- demand_and_capacity_by_zone$capacity\n
Run Code Online (Sandbox Code Playgroud)\n也就是说,我们不会移动已经在那里的人。
\n如果我们决定可以移动它们,规则就变成了Allocation >= 0
,我们就得到了埃尔文的答案。
coefs_con <- list() \nfor (z in zones) {\n coefs_con_zone <- list() \n for(categ in categories) {\n coefs_alloc <- alloc_coefs_template\n coefs_alloc[z, categ] <- 1\n # flatten the array\n coefs_con_zone[[categ]] <- c(move_coefs_template, coefs_alloc)\n }\n coefs_con[[z]] <- do.call(rbind, coefs_con_zone)\n}\n\n# stack the flattened arrays to build a matrix\nf.con3 <- do.call(rbind, coefs_con)\nf.dir3 <- rep(">=", n_zones * n_cat)\n\nif (move_new_only) {\n f.rhs3 <- c(t(population_mat))\n} else {\n f.rhs3 <- rep(0, n_zones * n_cat)\n}\n
Run Code Online (Sandbox Code Playgroud)\nf.con <- rbind(f.con1, f.con2, f.con3)\nf.dir <- c(f.dir1, f.dir2, f.dir3)\nf.rhs <- c(f.rhs1, f.rhs2, f.rhs3)\n
Run Code Online (Sandbox Code Playgroud)\n# compute the solution and visualize it in the array\nresults_raw <- lp("min", f.obj, f.con, f.dir, f.rhs)$solution\nresults_moves <- move_coefs_template\nresults_moves[] <- \n results_raw[1:length(results_moves)]\nresults_allocs <- alloc_coefs_template\nresults_allocs[] <- \n results_raw[length(results_moves)+(1:length(results_allocs))]\nresults_moves\n#> , , A\n#> \n#> 1 2 3 4\n#> 1 0 0 0 0\n#> 2 0 0 3 0\n#> 3 0 0 0 0\n#> 4 13 0 2 0\n#> \n#> , , B\n#> \n#> 1 2 3 4\n#> 1 0 0 0 0\n#> 2 0 0 0 0\n#> 3 0 0 0 0\n#> 4 0 0 57 0\n#> \n#> , , C\n#> \n#> 1 2 3 4\n#> 1 0 0 0 0\n#> 2 0 0 0 0\n#> 3 0 0 0 0\n#> 4 8 0 0 0\n#> \n#> , , D\n#> \n#> 1 2 3 4\n#> 1 0 0 0 0\n#> 2 0 0 0 0\n#> 3 0 0 0 0\n#> 4 0 38 0 0\n\nresults_allocs\n#> A B C D\n#> 1 151 99 168 47\n#> 2 142 83 81 87\n#> 3 139 178 33 34\n#> 4 76 139 58 58\n
Run Code Online (Sandbox Code Playgroud)\n# format as tidy data frame\nresults_df <-\n as.data.frame.table(results_moves) %>% \n setNames(c("from", "to", "category", "n")) %>% \n filter(n != 0) %>% \n mutate_at(c("from", "to"), as.numeric) %>% \n mutate_if(is.factor, as.character)\n\nresults_df\n#> from to category n\n#> 1 4 1 A 13\n#> 2 2 3 A 3\n#> 3 4 3 A 2\n#> 4 4 3 B 57\n#> 5 4 1 C 8\n#> 6 4 2 D 38\n
Run Code Online (Sandbox Code Playgroud)\npopulation_and_demand_by_category_and_zone <-\n bind_rows(\n results_df %>% \n group_by(category, zone = to) %>% \n summarize(correction = sum(n), .groups = "drop"),\n results_df %>% \n group_by(category, zone = from) %>% \n summarize(correction = -sum(n), .groups = "drop"),\n ) %>% \n left_join(population_and_demand_by_category_and_zone, ., by = c("category", "zone")) %>% \n replace_na(list(correction =0)) %>% \n mutate(new_population = demand + correction)\n\npopulation_and_demand_by_category_and_zone\n#> # A tibble: 16 \xc3\x97 6\n#> category zone population demand correction new_population\n#> <chr> <dbl> <dbl> <dbl> <dbl> <dbl>\n#> 1 A 1 115 138 13 151\n#> 2 A 2 121 145 -3.00 142\n#> 3 A 3 112 134 5.00 139\n#> 4 A 4 76 91 -15.0 76\n#> 5 B 1 70 99 0 99\n#> 6 B 2 59 83 0 83\n#> 7 B 3 86 121 57 178\n#> 8 B 4 139 196 -57 139\n#> 9 C 1 142 160 8 168\n#> 10 C 2 72 81 0 81\n#> 11 C 3 29 33 0 33\n#> 12 C 4 58 66 -8 58\n#> 13 D 1 22 47 0 47\n#> 14 D 2 23 49 38 87\n#> 15 D 3 16 34 0 34\n#> 16 D 4 45 96 -38 58\n\n\ndemand_and_capacity_by_zone <-\n population_and_demand_by_category_and_zone %>% \n group_by(zone) %>% \n summarise(population = sum(population), correction = sum(correction), new_population = sum(new_population)) %>% \n left_join(demand_and_capacity_by_zone, ., by = "zone")\n#> `summarise()` ungrouping output (override with `.groups` argument)\n \ndemand_and_capacity_by_zone\n#> # A tibble: 4 \xc3\x97 7\n#> zone demand capacity capacity_exceeded population correction new_population\n#> <dbl> <dbl> <dbl> <lgl> <dbl> <dbl> <dbl>\n#> 1 1 444 465 FALSE 349 21 465\n#> 2 2 358 393 FALSE 275 35 393\n#> 3 3 322 500 FALSE 243 62 384\n#> 4 4 449 331 TRUE 318 -118 331\n
Run Code Online (Sandbox Code Playgroud)\n我们看到人口从未减少并且一直处于容量不足的状态。
\n