\( \newcommand{\bm}[1]{\boldsymbol{\mathbf{#1}}} \)
library(magrittr)
library(readr)
coexi <- read_delim("data_coexi.txt", 
    "\t", escape_double = FALSE, col_names = FALSE, 
    trim_ws = TRUE) %>%
    set_names(c("time", "R1", "R2", "Predator"))

R1_win <- read_delim("data_R1_exclude_R2.txt", 
    "\t", escape_double = FALSE, col_names = FALSE, 
    trim_ws = TRUE) %>%
    set_names(c("time", "R1", "R2", "Predator"))

R2_win <- read_delim("data_R2_exclude_R1.txt", 
    "\t", escape_double = FALSE, col_names = FALSE, 
    trim_ws = TRUE) %>%
    set_names(c("time", "R1", "R2", "Predator"))

Model

\[\begin{equation} \begin{split} &\frac{d R_1}{dt} = [r_1 (1 - \frac{R_1}{K_1}) - a_1 P] ~ R_1 \\ &\frac{d R_2}{dt} = [r_2 (1 - \frac{R_2}{K_2}) - a_2 P] ~ R_2 \\ &\frac{d P}{dt} = (b_1 R_1 + b_2 R_2 - m) ~ P \end{split} \end{equation}\]

R1, R2 Coexist

double a1=0.7, a2=0.7, b1=1.0, b2=1.0;
double r1=0.8, r2=0.85, K1=1.2, K2=1.2, m=0.8;
library(ggplot2)
ggplot(coexi)+
    geom_point(aes(x=time, y=R1, color="R1"), size=0.01)+
    geom_point(aes(x=time, y=R2, color="R2"), size=0.01)+
    geom_point(aes(x=time, y=Predator, color="Predator"), size=0.01)+
    guides(colour = 
            guide_legend(override.aes =
                             list(size=2)))+
    scale_x_continuous(
        breaks = seq(0,300,50))+
    labs(x="Time", y="Abundance",color="",
         title=expression('R'[1]*' R'[2]*' Coexist'))