From ee5457b4a000278079962d9845f21ec27e5c5a10 Mon Sep 17 00:00:00 2001 From: Hannes Kuchelmeister Date: Fri, 23 Oct 2020 16:25:33 +0200 Subject: [PATCH] add r assignments from course --- assignment_1.r | 172 +++++++++++++++++++++++++++++++++++++ assignment_2.r | 126 +++++++++++++++++++++++++++ assignment_3.r | 227 +++++++++++++++++++++++++++++++++++++++++++++++++ 3 files changed, 525 insertions(+) create mode 100644 assignment_1.r create mode 100644 assignment_2.r create mode 100644 assignment_3.r diff --git a/assignment_1.r b/assignment_1.r new file mode 100644 index 0000000..0e660e3 --- /dev/null +++ b/assignment_1.r @@ -0,0 +1,172 @@ +library(DeliveryMan) + +create_node <- function(posX, posY, cost = .Machine$integer.max) { + return(list(x = posX, y = posY, cost = cost)) +} + +distance <- function (x1, y1, x2, y2){ + return(abs(x1 - x2) + abs(y1 - y2)) +} + +myFunction <- function(roads, car, packages) { + # roads-> list of two matrixes for traffic conditions hroads and vroads (coordintes bottom left to top right ) + # coordinates stat <1,1> and go to + # list of information about your car -> load -> mem -> nextMove + # matrix about packages + car$nextMove = 5 + + destX = 1 + destY = 1 + if(car$load == 0) { + # closest package + distances = distance(car$x, car$y, packages[,1], packages[,2]) + # set packages we already collected to infinite distance + distances[packages[, 5] != 0] = .Machine$integer.max + minDistanceIndex = which.min(distances) + destX = packages[minDistanceIndex, 1] + destY = packages[minDistanceIndex, 2] + + } else { + destX = packages[car$load, 3] + destY = packages[car$load, 4] + } + path = find_path(car$x, car$y, destX, destY, roads) + nextPoint = path[[1]] + car$nextMove = next_move(car$x, car$y, nextPoint$x, nextPoint$y) + return(car) +} + +get_neighbour <- function(posX, posY, cameFromX, cameFromY, goalX, goalY, costToNode) { + totalCost = costToNode + distance(posX, posY, goalX, goalY) + toVisit = TRUE + return(create_node_info_element(posX, posY, cameFromX, cameFromY, costToNode, totalCost, toVisit =TRUE)) +} + +get_neighbours <- function (posX, posY, goalX, goalY, costToParentNode, roads) { + maxY = length(roads$hroads[1,]) + maxX = length(roads$vroads[,1]) + + n1 = NULL + n2 = NULL + n3 = NULL + n4 = NULL + + if (posX > 1) { + cost = roads$hroads[posX-1, posY] + costToParentNode + n1 = create_node(posX-1,posY, cost) + } + if (posX < maxX) { + cost = roads$hroads[posX, posY] + costToParentNode + n2 = create_node(posX + 1, posY, cost) + } + if (posY > 1) { + cost = roads$vroads[posX, posY - 1] + costToParentNode + n3 = create_node(posX,posY - 1, cost) + } + if (posY < maxY) { + cost = roads$vroads[posX, posY] + costToParentNode + n4 = create_node(posX,posY + 1, cost) + } + return(list(n1, n2, n3, n4)) +} + +find_path <- function (startX, startY, goalX, goalY, roads) { + pq_cost = c() + pq_list = list() + + maxX = length(roads$vroads[,1]) + maxY = length(roads$hroads[1,]) + cameFromX <<- matrix(0, nrow = maxX, ncol = maxY) + cameFromY <<- matrix(0, nrow = maxX, ncol = maxY) + costToNode <<- matrix(.Machine$double.xmax, nrow = maxX, ncol = maxY) + visited <<- matrix(FALSE, nrow = maxX, ncol = maxY) + + maxY = length(roads$hroads[1,]) + maxX = length(roads$vroads[,1]) + + cost = 0 + node = create_node(startX, startY, 0) + # priority queue insert + insert.order <- findInterval(cost, pq_cost) + pq_cost <- append(pq_cost, cost, insert.order) + pq_list <- append(pq_list, list(node), insert.order) + + + costToNode[startX, startY] = 0 + cameFromX[startX, startY] = startX + cameFromY[startX, startY] = startY + + current = NULL + while (length(pq_cost) != 0) { + # pop queue + current <- pq_list[[1]] + pq_cost <- pq_cost[-1] + pq_list <- pq_list[-1] + + if (current$x == goalX & current$y == goalY) { + path_list = NULL + # walk backwards until starting element is reached + while (current$x != startX | current$y != startY) { + # get next element + path_list = append(list(current), path_list) + current = create_node(cameFromX[current$x, current$y], cameFromY[current$x, current$y], costToNode[current$x, current$y]) + } + return(path_list) + } + + neighbours = get_neighbours(current$x, current$y, goalX, goalY, costToNode[current$x, current$y], roads) + + for(n in neighbours) { + if(is.null(n) == FALSE) { + if(n$cost < costToNode[n$x, n$y]) { + cameFromX[n$x,n$y] = current$x + cameFromY[n$x,n$y] = current$y + + costToNode[n$x, n$y] = n$cost + + totalCost = n$cost + distance(n$x, n$y, goalX, goalY) + + # priority queue insert + insert.order <- findInterval(totalCost, pq_cost) + pq_cost <- append(pq_cost, totalCost, insert.order) + pq_list <- append(pq_list, list(n), insert.order) + } + } + } + + } + #Failed to find path + return(NULL) +} + +next_move <- function(carX, carY, dirX, dirY) { + if(!is.numeric(dirX) || !is.numeric(dirY)) { + return(0) + } + xDir = carX - dirX + yDir = carY - dirY + + if (xDir + yDir < 0) { + if (xDir == 0) { + # up + return(8) + } else { + # right + return(6) + } + } else if (xDir + yDir > 0){ + if (xDir == 0) { + # down + return(2) + } else { + # left + return(4) + } + } + return(0) +} + +startTime = Sys.time() +print(testDM(myFunction, timeLimit = 250)) +endTime = Sys.time() +print(endTime - startTime) diff --git a/assignment_2.r b/assignment_2.r new file mode 100644 index 0000000..1522748 --- /dev/null +++ b/assignment_2.r @@ -0,0 +1,126 @@ +library(WheresCroc) + +myFunction = function (moveInfo, readings, positions, edges, probs){ + checkIfInTop = 3 + # reset state if new game + if((moveInfo$mem$status >= 0)) { # status new game has started + #reset states + dim = length(probs$salinity[,1]) + + moveInfo$mem$state = rep(1/dim, dim) + if ((moveInfo$mem$status == 0)) { # status first game has started + connections = diag(1, nrow = dim, ncol = dim) + + # set connections in matrix + connections[edges] = 1 + connections[edges[, c(2,1)]] = 1 + + # calculate transition_matrix based on connections + divide_by = colSums(connections) + moveInfo$mem$connections = connections + moveInfo$mem$transition_matrix = sweep(connections, 2, divide_by, FUN = '/') + moveInfo$mem$checkIfInTop = 22 + } + + moveInfo$mem$status = -1 + } + + connections = moveInfo$mem$connections + state = moveInfo$mem$state + transition_matrix = moveInfo$mem$transition_matrix + player_pos = positions[3] + + # update sate based on backpackers + if (is.na(positions[1]) == FALSE) { + if (positions[1] < 0) { # gets eaten + state[] = 0 + state[positions[1]] = 1 + } else { + state[positions[1]] = 0 + } + } + if (is.na(positions[2]) == FALSE) { + if (positions[2] < 0) { # gets eaten + state[] = 0 + state[positions[2]] = 1 + } else { + state[positions[2]] = 0 + } + } + + probs_salinity = dnorm(readings[1], mean = probs$salinity[,1], sd = probs$salinity[,2]) + probs_phosphate = dnorm(readings[2], mean = probs$phosphate[,1], sd = probs$phosphate[,2]) + probs_nitrogen = dnorm(readings[3], mean = probs$nitrogen[,1], sd = probs$nitrogen[,2]) + + # calculating the emission matrix + emission_matrix = diag(probs_salinity * probs_phosphate * probs_nitrogen) + + # probabilityies of crocs position next turn + state = emission_matrix %*% state + + # normalise state + state = state / sum(state) + + node_to_move_to = which.max(state) # most likely next round + path = c(0) + if (node_to_move_to != player_pos) { + path = findPath(player_pos, node_to_move_to, connections) + } else { + order = order(state, decreasing = TRUE) + node_to_move_to_plan_b = order[2] + path = append(node_to_move_to, findPath(player_pos, node_to_move_to_plan_b, connections)) + } + ranks = rank(state) + if((ranks[player_pos] > (length(ranks) - checkIfInTop)) & (ranks[path[1]] < ranks[player_pos] )) { + moveInfo$moves = c(0, path[1]) + state[player_pos] = 0 + } else if(ranks[path[1]] > (length(ranks) - checkIfInTop)){ + moveInfo$moves = c(path[1], 0) + state[path[1]] = 0 + } else { + moveInfo$moves = c(path[1], path[2]) + } + + state = transition_matrix %*% state + + moveInfo$mem$state = state + return(moveInfo) +} + +findPath <- function(start, goal, matrix) { + nodesNumber = length(matrix[1,]) + q = c(start) + visited = rep(F, nodesNumber) + visited[start] = T + prev = rep(NULL, nodesNumber) + + while(length(q) > 0) { + current = q[1] + q = q[-1] + + if (current == goal) { + path = c(goal) + node = prev[goal] + while (!is.na(prev[node])) { + path = c(node, path) + node = prev[node] + } + return(path) + } + + neighbours = which(matrix[current,] > 0) + + for (neighbour in neighbours) { + if (!visited[neighbour]) { + q = c(q, neighbour) + visited[neighbour] = T + prev[neighbour] = current + } + } + } + return(NA) +} + +#runWheresCroc(myFunction, showCroc = TRUE) +testWC(myFunction, verbose = 1) +#runWheresCroc(WheresCroc::manualWC) diff --git a/assignment_3.r b/assignment_3.r new file mode 100644 index 0000000..b7c5ee8 --- /dev/null +++ b/assignment_3.r @@ -0,0 +1,227 @@ +library(Diagnostics) + + + +learn <- function(cases) { + # order: P(Pn), P(Te |Pn), P(VTB), P(TB| VTB), P(Sm), P(LC |Sm), P(BR |Sm), P(XR |Pn, VTB, LC), P(Dy |LC, BR) + pn = table(cases$Pn) + 1 + pn = pn / sum(pn) + + #special handling for temperature + te__pn = generate_te_model(cases) + + vtb = table(cases$VTB) + 1 + vtb = vtb / sum(vtb) + + tb__vtb = table(cases$VTB, cases$TB) + 1 + for (i in 1:2) { + tb__vtb[i, ] = tb__vtb[i, ] / sum(tb__vtb[i, ]) + } + + sm = table(cases$Sm) + 1 + sm = sm / sum(sm) + + lc__sm = table(cases$Sm, cases$LC) + 1 + for (i in 1:2) { + lc__sm[i, ] = lc__sm[i, ] / sum(lc__sm[i, ]) + } + + br__sm = table(cases$Sm, cases$Br) + 1 + for (i in 1:2) { + br__sm[i, ] = br__sm[i, ] / sum(br__sm[i, ]) + } + + xr__pn_tb_lc = table(cases$Pn, cases$TB, cases$LC, cases$XR) + 1 + for (i in 1:2) { + for (j in 1:2) { + for (k in 1:2) { + xr__pn_tb_lc[i, j, k, ] = xr__pn_tb_lc[i, j, k, ] / sum(xr__pn_tb_lc[i, j, k, ]) + } + } + } + + dy__lc_br = table(cases$LC, cases$Br, cases$Dy) + 1 + for (i in 1:2) { + for (j in 1:2) { + dy__lc_br[i, j, ] = dy__lc_br[i, j, ] / sum(dy__lc_br[i, j, ]) + } + } + + model = list() + model$pn = pn + model$te__pn = te__pn + model$vtb = vtb + model$tb__vtb = tb__vtb + model$sm = sm + model$lc__sm = lc__sm + model$br__sm = br__sm + model$xr__pn_tb_lc = xr__pn_tb_lc + model$dy__lc_br = dy__lc_br + + return(model) +} + +generate_te_model <- function(cases) { + cases0 = list() + cases1 = list() + + for (i in 1:length(cases$Te)) { + if (cases$Pn[i] == 0) { + cases0 = append(cases0, cases$Te[i]) + } else{ + cases1 = append(cases1, cases$Te[i]) + } + } + cases0 = as.vector(unlist(cases0)) + cases1 = as.vector(unlist(cases1)) + + l0 = list() + l0$mean = mean(cases0) + l0$sd = sd(cases0) + l1 = list() + l1$mean = mean(cases1) + l1$sd = sd(cases1) + + return(list(l0, l1)) # access with var_name[[1]]$sd, var_name[[2]]$sd, etc +} + +diagnose <- function(model, cases) { + samples_number = 1500 + burn_period = 0.1 + + for (i in 1:length(cases[,1])){ + case = cases[i,] + samples = generate_samples(cases[i,], samples_number, burn_period, model) + case = make_predictions(case, samples) + cases[i,] = case + } + return((c(cases$Pn, cases$TB, cases$LC, cases$Br))) +} + +evaluate_var_prob <- function(samples, var_name) { + samples_count = length(samples[,1]) + var_column = samples[var_name] + var_true_count = sum(var_column) + true_prob = var_true_count / samples_count + return(true_prob) +} + +make_predictions <- function(case, samples) { + + case["Pn"] = evaluate_var_prob(samples, "Pn") + case["TB"] = evaluate_var_prob(samples, "TB") + case["LC"] = evaluate_var_prob(samples, "LC") + case["Br"] = evaluate_var_prob(samples, "Br") + + return(case) +} + +invert <- function(val) { + if (val == 0) { + return(1) + } else { + return(0) + } +} + +calc_conditional_probs <- function(sample, model) { + Pn = sample$Pn + Te = sample$Te + VTB = sample$VTB + TB = sample$TB + Sm = sample$Sm + LC = sample$LC + Br = sample$Br + XR = sample$XR + Dy = sample$Dy + return( + model$pn[Pn + 1] * + dnorm( + mean = model$te__pn[[Pn + 1]]$mean, + sd = model$te__pn[[Pn + 1]]$sd, + x = Te + ) * + model$vtb[VTB + 1] * + model$tb__vtb[VTB + 1, TB + 1] * + model$sm[Sm + 1] * + model$lc__sm[Sm + 1, LC + 1] * + model$br__sm[Sm + 1, Br + 1] * + model$xr__pn_tb_lc[Pn + 1, TB + 1, LC + 1, XR + 1] * + model$dy__lc_br[LC + 1, Br + 1, Dy + 1] + ) +} + +propose_value_for_unknown <- function(sample, unknown_var, model) { + old_sample = sample + new_sample = sample + new_sample[unknown_var] = invert(new_sample[unknown_var]) + p_old = calc_conditional_probs(old_sample, model) + p_new = calc_conditional_probs(new_sample, model) + + if (p_new > p_old) { + return(new_sample) + } else { + threshold = p_new / p_old + random = runif(1, 0, 1) + if (random < threshold) { + return(new_sample) + } else { + return(old_sample) + } + } + +} + +generate_one_sample <- function(samples, number, model) { + sample = samples[number,] + if (number == 1) { #first sample, assign unknown variables to random values + sample$Pn = sample(0:1, 1) + sample$TB = sample(0:1, 1) + sample$LC = sample(0:1, 1) + sample$Br = sample(0:1, 1) + } else { + prev_sample = samples[number - 1,] + sample$Pn = prev_sample$Pn + sample$TB = prev_sample$TB + sample$LC = prev_sample$LC + sample$Br = prev_sample$Br + } + + sample = propose_value_for_unknown(sample, "Pn", model) + sample = propose_value_for_unknown(sample, "TB", model) + sample = propose_value_for_unknown(sample, "LC", model) + sample = propose_value_for_unknown(sample, "Br", model) + + return(sample) +} + +generate_samples <- function(case, samples_number, burn_period, model) { + + samples = hist + samples = samples[1:samples_number,] # samples_number cannot be more than hist size, maybe fix it later + samples$Te = case$Te + samples$VTB = case$VTB + samples$Sm = case$Sm + samples$XR = case$XR + samples$Dy = case$Dy + + for (i in 1:samples_number) { + samples[i,] = generate_one_sample(samples, i, model) + } + + samples_burned = samples_number * burn_period + samples = samples[samples_burned:length(samples[,1]),] + + return(samples) +} +#run_count = 10 +#runs = vector(length = run_count) +runDiagnostics(learn, diagnose, verbose = 2) +#for (i in 1:run_count) { +# runs[i] = runDiagnostics(learn, diagnose, verbose = 2) +#} +# print("Result: ") +# cat("mean: ", mean(runs), "\n") +# cat("sd: ", sd(runs), "\n") +# cat("worse than 0.01345: ", length(runs[runs > 0.01345]), "\n") +# cat("worse than 0.013375: ", length(runs[runs > 0.013375]), "\n") \ No newline at end of file