per_variable_censoring.jl

# Install the feature branch of Copulas.jl : 
# ] add https://github.com/lrnv/Copulas.jl#lrnv/issue345

using Copulas, Distributions

θ  = 2.0
C  = ClaytonCopula(2, θ)
X₁ = Exponential(1.0)          # margin 1 — observed
X₂ = Exponential(1.0)          # margin 2 — right-censored at c2
x₁ = 0.7                       # observed event time (dim 1)
c₂ = 1.3                       # dim 2 right-censored: we only know T2 > c2

u₁, u₂ = cdf(X₁, x₁), cdf(X₂, c₂)

X = SklarDist(C, (X₁, X₂))

# You want to compute logP(X₁ = x₁, X₂ ≤ c₂), which is indeed your log-likelihood. 
# Using the same decomposition as you with the conditional : 
# logP(X₁ = x₁, X₂ ≤ c₂) = logP(X₂ ≤ c₂ | X₁ = x₁) + logP(X₁ = x₁)

# You can literally implement this decomposition as :
logcdf(condition(X, 1, x₁), c₂) + logpdf(X₁, x₁) # = -1.0049...

# Or equivalently if you move the conditional to the copula scale: 
logcdf(condition(C, 1, u₁), u₂) + logpdf(X₁, x₁) # = -1.0049...

# while for non-censoreed observations you can simply use : 
logpdf(X, [x₁, c₂])

# A more generic version of per-variable censored likelyhood might look like: 

function _split(x, δ)
    d, n = length(x), sum(δ)
    i = findall(Bool.(δ))
    xᵢ  = n == 1   ? x[i[1]]   : x[δ]
    x₋ᵢ = n == d-1 ? x[.!δ][1] : x[.!δ]
    return i, xᵢ, x₋ᵢ
end
function per_var_censored_llh(X, x, δ)

    all(δ .== 0) && return logcdf(X, x) # everything is censored
    all(δ .== 1) && return logpdf(X, x) # everything is observed
    i, xᵢ, x₋ᵢ = _split(x, δ)
    return logcdf(condition(X, i, xᵢ), x₋ᵢ) + logpdf(subsetdims(X, i), xᵢ)
end

# Let's try it: 
X = SklarDist(ClaytonCopula(7, 0.3), (X₁, X₁, X₁, Normal(), X₂, X₂, X₂))
x = rand(X)
δ = Bool.(rand(Bernoulli(0.8), 7))
per_var_censored_llh(X, x, δ) # Looks correct to me ! 

# And this shoul work with any model supporting our interface. 
# I'm not sure about performances tho !