I am working with the R programming language.

- Suppose there are 100 buses (bus_1, bus_2,... bus_100).
- The bus company sends the first 5 buses (bus_1, bus_2... bus_5) on the first day.
- On the first day, each bus has a 0.5 probability of breaking down (once a bus breaks down, it is out of service permanently).
- If the bus does not break down on the first day, its sent out on the second day but now it has a probability of breaking down of 0.49. That is each day a bus doesn't break down, the probability of breaking down reduces by 0.01.
- When a bus breaks down, it is replaced with the next bus (e..g if bus1, bus2, bus3,bus4, bus5 are there bus2 and bus4 break on the same day, the next day we will have bus1, bus6, bus3, bus7 , bus5).
- At any given day, there can only be a max of 5 buses out on the road ... towards the end, there might be 4 buses, 3 buses ... until they all break down.

**My Question:** I am trying to write a simulation which simulates this situation until all buses break down? The final result should be a data frame with 11 columns: day_number (1,2,3...n), x1, x2, x3, x4, x5 (these represent the bus number in that position), p1,p2,p3,p4,p5 (these represent the probabilities of each bus breaking down on that day's row).

Here is what I tried so far using basic loops:

```
bus_count <- 100
bus_prob <- rep(0.5, bus_count)
bus_active <- 1:5
results <- data.frame(matrix(ncol = 11, nrow = 0))
colnames(results) <- c("day_number", "x1", "x2", "x3", "x4", "x5", "p1", "p2", "p3", "p4", "p5")
day_number <- 0
while(length(bus_active) > 0) {
day_number <- day_number + 1
breakdown <- runif(length(bus_active)) < bus_prob[bus_active]
bus_prob[bus_active] <- pmax(bus_prob[bus_active] - 0.01, 0)
next_bus <- (max(bus_active)+1):(max(bus_active)+sum(breakdown))
next_bus <- next_bus[next_bus <= bus_count]
bus_prob[next_bus] <- 0.5
bus_active <- c(bus_active[!breakdown], next_bus)
results <- rbind(results, c(day_number, bus_active, bus_prob[bus_active]))
}
# rename bus values to actual names
results[,2:6] <- apply(results[,2:6], 2, function(x) paste0("bus_", x))
print(results)
```

This is not correct. I am noticing the following errors (column names did not come out correct):

Problem 1: The bus ordering is incorrect. On the first day, the buses should be bus1, bus2, bus3, bus4, bus5.

```
X1 X3 X5 X6 X7 X8 X0.49 X0.49.1 X0.5 X0.5.1 X0.5.2
1 1 bus_3 bus_5 bus_6 bus_7 bus_8 0.49 0.49 0.50 0.50 0.5
2 2 bus_3 bus_6 bus_9 bus_10 bus_11 0.48 0.49 0.50 0.50 0.5
```

Problem 2: Buses are "coming back from the dead". E.g. Day 15 vs Day 47, bus47 comes back

```
15 15 bus_38 bus_45 bus_46 bus_47 bus_48 0.46 0.49 0.50 0.50 0.5
47 47 bus_47 bus_47 bus_47 bus_47 bus_47 47.00 47.00 47.00 47.00 47.0
```

Problem 3: The same bus appears multiple times in the same row:

```
47 47 bus_47 bus_47 bus_47 bus_47 bus_47 47.00 47.00 47.00 47.00 47.0
```

Problem 4: The probabilities are not in the correct ranges (e.g. can only be between 0.5 and 0)

```
37 37 bus_98 bus_99 bus_98 bus_0.5 bus_0.5 0.50 37.00 98.00 99.00 98.0
38 38 bus_99 bus_98 bus_100 bus_0.49 bus_0.49 0.50 38.00 99.00 98.00 100.0
```

Can someone please help me fix these problems and write this code correctly?

Thanks!