Yes, indexing the model output and target should work and the gradient would be backpropagated to the selected values:
output = torch.randn(10, 1, requires_grad=True)
target = torch.randint(0, 2, (10, 1)).float()
criterion = nn.BCEWithLogitsLoss()

batch_idx = torch.tensor([1, 3, 5, 7])
loss = criterion(output[batch_idx], target[batch_idx])
loss.backward()

print(output.grad)
> tensor([[ 0.0000],
          [-0.1686],
          [ 0.0000],
          [ 0.1138],
          [ 0.0000],
          [ 0.2076],
          [ 0.0000],
          [-0.1147],
          [ 0.0000],
          [ 0.0000]])