Created
October 2, 2024 07:21
-
-
Save lesserfish/eaef33757ee4bf6c19096d4e723f917b to your computer and use it in GitHub Desktop.
Fast Michaelis-Menten equations
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| # Make sure to have Rato installed on your system! | |
| # Michaelis Menten equation | |
| michaelis.menten <- function(t, x, parms) { | |
| f <- parms$f | |
| h <- parms$h | |
| B <- parms$B | |
| beta <- parms$beta | |
| dxdt <- -B*(x^f) + beta*(x^h)/(x^h + 1) # Exemplo de equação: dy/dt = -0.1 * y | |
| list(dxdt) | |
| } | |
| x_eff <- function(M, x) { | |
| indegree <- colSums(M) | |
| outdegree <- rowSums(M) | |
| n <- nrow(M) | |
| ones <- matrix(1, nrow = n, ncol = 1) | |
| x_eff <- (outdegree * x) / sum(indegree) | |
| num <- as.numeric(t(ones) %*% M %*% x) | |
| den <- as.numeric(t(ones) %*% M %*% ones) | |
| return(num / den) | |
| } | |
| fetch_parameters <- function(M){ | |
| indegree <- colSums(M) | |
| outdegree <- rowSums(M) | |
| initial_params <- list('f' = 1, 'h' = 2, 'B' = 0.1) | |
| initial_params$beta <- ((indegree %*% outdegree) / sum(indegree)) | |
| return(initial_params) | |
| } | |
| deform <- function(M, level){ | |
| index <- sample(1:nrow(M), level) | |
| M[index, ] <- 0 | |
| M[, index] <- 0 | |
| return(M) | |
| } | |
| lerpcol <- function(a, b, l){ | |
| o <- l * b + (1 - l) * a | |
| print(o) | |
| return(rgb(o[1], o[2], o[3])) | |
| } | |
| # Load the graph! | |
| g <- Rato::graph.from.csv("./Demos/IL17/IL17.nodes.csv", "./Demos/IL17/IL17.edges.csv", sep=",", header=TRUE) | |
| # Initial plot | |
| plot(NA, NA, , xlab="Time", ylab="Overall gene expression", xlim=c(0, 100), ylim = c(0, 65), main="Fast Michaelis Menten equations") | |
| for(level in 1:nrow(g$M)){ | |
| col = lerpcol(c(1, 0, 0), c(0, 0, 1), level/nrow(g$M)) | |
| # Perturbation function | |
| M <- deform(g$M, level) | |
| # Calculate the initial value | |
| xeff <- x_eff(M, g$initial_values) | |
| # Solve the differential equation | |
| solution <- deSolve::ode(y = xeff, func = michaelis.menten, parms = fetch_parameters(M), times = seq(0, 100, 0.1)) | |
| # Plot the trajectory | |
| lines(x = solution[, 1], y = solution[, 2], xlab="Time", ylab="Overall gene expression", type="l", col=col, lwd=2) | |
| } | |
| legend("bottomright", # Position of the legend | |
| legend = c("No perturbation", "Fully perturbed"), # Legend labels | |
| col = c("red", "blue"), # Colors for lines | |
| lty = 1, # Line type (solid) | |
| lwd = 2, # Line width | |
| ) | |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment