def bgplvm_test_model(optimize=False, verbose=1, plot=False, output_dim=200, nan=False):
    import GPy

    num_inputs = 13
    num_inducing = 5
    if plot:
        output_dim = 1
        input_dim = 3
    else:
        input_dim = 2
        output_dim = output_dim

    # generate GPLVM-like data
    X = _np.random.rand(num_inputs, input_dim)
    lengthscales = _np.random.rand(input_dim)
    k = GPy.kern.RBF(input_dim, .5, lengthscales, ARD=True)
    K = k.K(X)
    Y = _np.random.multivariate_normal(_np.zeros(num_inputs), K, (output_dim,)).T

    # k = GPy.kern.RBF_inv(input_dim, .5, _np.ones(input_dim) * 2., ARD=True) + GPy.kern.bias(input_dim) + GPy.kern.white(input_dim)
    # k = GPy.kern.linear(input_dim)# + GPy.kern.bias(input_dim) + GPy.kern.white(input_dim, 0.00001)
    # k = GPy.kern.RBF(input_dim, ARD = False)  + GPy.kern.white(input_dim, 0.00001)
    # k = GPy.kern.RBF(input_dim, .5, _np.ones(input_dim) * 2., ARD=True) + GPy.kern.RBF(input_dim, .3, _np.ones(input_dim) * .2, ARD=True)
    # k = GPy.kern.RBF(input_dim, .5, 2., ARD=0) + GPy.kern.RBF(input_dim, .3, .2, ARD=0)
    # k = GPy.kern.RBF(input_dim, .5, _np.ones(input_dim) * 2., ARD=True) + GPy.kern.linear(input_dim, _np.ones(input_dim) * .2, ARD=True)

    p = .3

    m = GPy.models.BayesianGPLVM(Y, input_dim, kernel=k, num_inducing=num_inducing)

    if nan:
        m.inference_method = GPy.inference.latent_function_inference.var_dtc.VarDTCMissingData()
        m.Y[_np.random.binomial(1, p, size=(Y.shape)).astype(bool)] = _np.nan
        m.parameters_changed()

    #===========================================================================
    # randomly obstruct data with percentage p
    #===========================================================================
    # m2 = GPy.models.BayesianGPLVMWithMissingData(Y_obstruct, input_dim, kernel=k, num_inducing=num_inducing)
    # m.lengthscales = lengthscales

    if plot:
        import matplotlib.pyplot as pb
        m.plot()
        pb.title('PCA initialisation')
        # m2.plot()
        # pb.title('PCA initialisation')

    if optimize:
        m.optimize('scg', messages=verbose)
        # m2.optimize('scg', messages=verbose)
        if plot:
            m.plot()
            pb.title('After optimisation')
            # m2.plot()
            # pb.title('After optimisation')

    return m
model for testing purposes. Samples from a GP with rbf kernel and learns
dataset/ETHPy150Open SheffieldML/GPy/GPy/examples/dimensionality_reduction.py/bgplvm_test_model