Non homogeneous Semi Markov data simulation for continuous time

[purplew2]

I need to simulate data for a non homogeneous Semi Markov process for continuous time. I've tried some codes in R but they seem to not be working properly. This is the code that I tried in R. Does anyone have any suggestions of making it better? Also, is there any packages for this specific task (Python or R)? I know SMM is for descrete time.


# Define the state labels
state_labels <- c("A", "B", "C", "D")

# Define the transition matrix with dimnames
trans_mat <- matrix(c(0, 0.5, 0.25, 0.25,
0, 0, 0.5, 0.5,
0, 0, 0, 1,
0.1, 0, 0, 0.9), nrow = 4, byrow = TRUE,
dimnames = list(state_labels, state_labels))

# Set the initial state
initial_state <- "A"

# Set the simulation time

# Function to simulate semi-Markov process
simulate_semi_markov <- function(trans_mat, state_labels, initial_state, sim_time) {
n_states <- length(state_labels)
time <- 0
state <- initial_state
states <- c(state)

while (time < sim_time) {
# Simulate next state transition
trans_probs <- trans_mat[state, ]
trans_probs[state] <- 0
trans_probs[trans_probs <= 0] <- 0
trans_probs <- trans_probs / sum(trans_probs)

next_state <- sample(state_labels, size = 1, prob = trans_probs)
rate <- trans_mat[state, next_state]

if (rate == 0) {
dwell_time <- sim_time - time
} else {
dwell_time <- rexp(1, rate = 1/rate)
if (is.na(dwell_time)) {
break
}
if (time + dwell_time > sim_time) {
dwell_time <- sim_time - time
}
}

time <- time + dwell_time
state <- next_state
states <- c(states, state)
}

return(list(time = time, states = states))
}

# Simulate the semi-Markov process
sim_data <- simulate_semi_markov(trans_mat, state_labels, initial_state, sim_time)

# Print the simulated data
print(sim_data)

# Extract the simulated states
simulated_states <- sim_data$states

# Print the simulated states
print(simulated_states)

 

[BigBeachBanana]

There are no specific packages in R or Python that provide built-in functions for simulating non-homogeneous semi-Markov processes. However, you can use general-purpose packages such as markovchain in R or Scipy in Python to build custom simulations.

Try this. Made some edits to your code to improve performance.

Spoiler

# Define the state labels
state_labels <- c("A", "B", "C", "D")

# Define the transition matrix with dimnames
trans_mat <- matrix(c(0, 0.5, 0.25, 0.25,
                      0, 0, 0.5, 0.5,
                      0, 0, 0, 1,
                      0.1, 0, 0, 0.9),
                    nrow = 4, byrow = TRUE,
                    dimnames = list(state_labels, state_labels))

# Set the initial state
initial_state <- "A"

# Set the simulation time
sim_time <- 10

# Function to simulate semi-Markov process
simulate_semi_markov <- function(trans_mat, state_labels, initial_state, sim_time) {
  n_states <- length(state_labels)
  time <- 0
  state <- initial_state
  states <- c(state)
 
  while (time < sim_time) {
    # Simulate next state transition
    trans_probs <- trans_mat[state, ]
    trans_probs[state] <- 0
    trans_probs <- trans_probs[trans_probs > 0]
    trans_probs <- trans_probs / sum(trans_probs)
  
    next_state <- sample(state_labels, size = 1, prob = trans_probs, replace = TRUE)
    rate <- trans_mat[state, next_state]
  
    if (rate == 0) {
      dwell_time <- sim_time - time
    } else {
      dwell_time <- rexp(1, rate = 1/rate)
      if (is.na(dwell_time)) {
        break
      }
      if (time + dwell_time > sim_time) {
        dwell_time <- sim_time - time
      }
    }
  
    time <- time + dwell_time
    state <- next_state
    states <- c(states, state)
  }
 
  return(list(time = time, states = states))
}

# Simulate the semi-Markov process
sim_data <- simulate_semi_markov