Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Error in quantile.default #26

Closed
gt7901b opened this issue Jun 22, 2021 · 3 comments
Closed

Error in quantile.default #26

gt7901b opened this issue Jun 22, 2021 · 3 comments

Comments

@gt7901b
Copy link

gt7901b commented Jun 22, 2021

Hi:

Thanks for developing this nice package.

I tried this command but got the

dsb_norm_prot = DSBNormalizeProtein(
cell_protein_matrix = positive_adt_matrix,
empty_drop_matrix = neg_adt_matrix,
denoise.counts = TRUE,
use.isotype.control = TRUE,
isotype.control.name.vec = rownames(citeseq_counts)[162:165])

Error in quantile.default(x, seq(from = 0, to = 1, length = n)): missing values and NaN's not allowed if 'na.rm' is FALSE
Traceback:

  1. DSBNormalizeProtein(cell_protein_matrix = positive_adt_matrix,
    . empty_drop_matrix = neg_adt_matrix, denoise.counts = TRUE,
    . use.isotype.control = TRUE, isotype.control.name.vec = rownames(citeseq_counts)[162:165])
  2. apply(norm_adt, 2, function(x) {
    . g = mclust::Mclust(x, G = 2, warn = F, verbose = F)
    . return(g$parameters$mean[1])
    . })
  3. FUN(newX[, i], ...)
  4. mclust::Mclust(x, G = 2, warn = F, verbose = F)
  5. eval(mc, parent.frame())
  6. eval(mc, parent.frame())
  7. mclustBIC(data = structure(c(-0.756590760281651, 0.706345357143731,
    . -0.571362010259017, 0.455675113554685, 0.178004789840622, 0.11274129347221,
    ...
    . NaN, NaN, -0.20047762928513, NaN, -0.0629144943331455, NaN, -0.0719656170716554,
    . -0.104914014711862, NaN), .Dim = c(217L, 1L), .Dimnames = list(
    . c("anti-human-CD140a-PDGFRalpha-TotalSeqC", "anti-human-CD140b-PDGFRbeta-TotalSeqC",
    . "anti-human-CD119-IFN-gamma-R-alpha-chain-TotalSeqC", "Rat-IgG1-kapa-isotype-Ctrl-TotalSeqC",

What do you think might be the reason?

before this command, I ran
Idents(seurat_obj) = "HTO_classification.global"

subset empty drop/background and cells

neg_object = subset(seurat_obj, idents = "Negative")
singlet_object = subset(seurat_obj, idents = "Singlet")
neg_adt_matrix = GetAssayData(neg_object, assay = "ADT", slot = 'counts') %>% as.matrix()
positive_adt_matrix = GetAssayData(singlet_object, assay = "ADT", slot = 'counts') %>% as.matrix()

the out put of rownames(citeseq_counts)[162:165] is
'Rat-IgG2b-kapa-Isotype-Ctrl-TotalSeqC''Mouse-IgG2b-kapa-isotype-Ctrl-TotalSeqC''Mouse-IgG2a-kapa-isotype-Ctrl-TotalSeqC''Mouse-IgG1-kapa-isotype-Ctrl-TotalSeqC'

thanks for your time

@MattPM
Copy link
Collaborator

MattPM commented Jun 23, 2021

Hi, please see: #6

This error occurs during denoising, (denoise = TRUE) when you have antibodies with 0 counts or close to 0 across all cells. To get rid of this error, check the distributions of the antibodies with e.g:
apply(positive_adt_matrix, 1, quantile)
apply(positive_adt_matrix, 1, max)

to find the protein(s) with basically no counts, then remove these from BOTH the background drops neg_adt_matrix and the cells positive_adt_matrix

Also see here for some more code examples:
https://github.com/niaid/dsb#step3 --> scroll down to "Optional step; remove proteins without staining"

Overall with the larger panels, it is better to remove those proteins with basically no data up front, it will just add noise to keep them in downstream analysis.

@gt7901b that should work, let me know if that fixes it.

@gt7901b
Copy link
Author

gt7901b commented Jun 23, 2021

Thank you so much. I was able to remove the proteins with no counts by doing

d1 = data.frame(pmax = apply(positive_adt_matrix, 1, max)) %>%
rownames_to_column('prot') %>% arrange(pmax) %>% head()
prot_names = rownames(positive_adt_matrix)
positive_adt_matrix = positive_adt_matrix[!prot_names == 'CD34-TotalSeqC', ]

No more error!

Another question regarding the isotype control
Do we need to know isotype control is for which anitbody?

thanks for your time

@MattPM
Copy link
Collaborator

MattPM commented Jun 24, 2021

Great, glad it worked. You do not need to know which antibody matches which isotype control. The isotypes are combined with the modeled background for the cell into the cell's technical component, which is then regressed out of the counts.

@MattPM MattPM closed this as completed Aug 4, 2021
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants