We prove the existence of positive solutions for the system
$$\begin{align*}\begin{cases}-\Delta_{p} u =\lambda a(x){f(v)}{u^{-\alpha}},\qquad x\in \Omega,\\-\Delta_{q} v = \lambda b(x){g(u)}{v^{-\beta}},\qquad x\in \Omega,\\u = v =0, \qquad x\in\partial \Omega,\end{cases}\end{align*}$$ where $\Delta_{r}z={\rm div}(|\nabla z|^{r-2}\nabla z)$, for $r>1$ denotes the r-Laplacian operator and $\lambda$ is a positive parameter, $\Omega$ is a bounded domain in $\mathbb{R}^{n}$, $n\geq1$ with sufficiently smooth boundary and $\alpha, \beta \in (0,1).$ Here $ a(x)$ and $ b(x)$ are $C^{1}$ sign-changing
functions that maybe negative near the boundary and $f,g $ are $C^{1}$ nondecreasing functions, such that $f, g :\ [0,\infty)\to [0,\infty);$ $f(s)>0,$ $g(s)>0$ for $s> 0$, $\lim_{s\to\infty}g(s)=\infty$ and
$$\lim_{s\to\infty}\frac{f(Mg(s)^{\frac{1}{q-1}})}{s^{p-1+\alpha}}=0,\qquad \forall M>0.$$
We discuss the existence of positive weak solutions when $f$, $g$, $a(x)$ and $b(x)$ satisfy certain additional conditions. We employ the method of sub-supersolution to obtain our results.