Modeling high-dimensional multivariate distributions is a computationally challenging task. Bayesian networks have been successfully used to reduce the complexity and simplify the problem with discrete variables. However, it lacks of a general model for continuous variables. In order to overcome this problem, Elidan (2010) proposed the model of copula bayesian networks (CBN) that reparametrizes bayesian networks with conditional copula functions. We propose a new learning algorithm for CBN based on a PC algorithm and a conditional independence test proposed by Bouezmarni, Rombouts, Taamouti (2009). This test being non-parametric, no model assumptions are made allowing it to be as general as possible. This algorithm is compared on generated data with the score based method proposed by Elidan (2010)}. Not only it proves to be faster, but also it generalizes well on data generated from distributions far from the gaussian model.