diff --git a/.Rbuildignore b/.Rbuildignore new file mode 100644 index 0000000..91114bf --- /dev/null +++ b/.Rbuildignore @@ -0,0 +1,2 @@ +^.*\.Rproj$ +^\.Rproj\.user$ diff --git a/.github/.gitignore b/.github/.gitignore new file mode 100644 index 0000000..2d19fc7 --- /dev/null +++ b/.github/.gitignore @@ -0,0 +1 @@ +*.html diff --git a/.github/workflows/R-CMD-check.yaml b/.github/workflows/R-CMD-check.yaml new file mode 100644 index 0000000..e36deaa --- /dev/null +++ b/.github/workflows/R-CMD-check.yaml @@ -0,0 +1,278 @@ +## Read more about GitHub actions the features of this GitHub Actions workflow +## at https://lcolladotor.github.io/biocthis/articles/biocthis.html#use_bioc_github_action +## +## For more details, check the biocthis developer notes vignette at +## https://lcolladotor.github.io/biocthis/articles/biocthis_dev_notes.html +## +## You can add this workflow to other packages using: +## > biocthis::use_bioc_github_action() +## +## Using GitHub Actions exposes you to many details about how R packages are +## compiled and installed in several operating system.s +### If you need help, please follow the steps listed at +## https://github.com/r-lib/actions#where-to-find-help +## +## If you found an issue specific to biocthis's GHA workflow, please report it +## with the information that will make it easier for others to help you. +## Thank you! + +## Acronyms: +## * GHA: GitHub Action +## * OS: operating system + +on: + push: + pull_request: + +name: R-CMD-check-bioc + +## These environment variables control whether to run GHA code later on that is +## specific to testthat, covr, and pkgdown. +## +## If you need to clear the cache of packages, update the number inside +## cache-version as discussed at https://github.com/r-lib/actions/issues/86. +## Note that you can always run a GHA test without the cache by using the word +## "/nocache" in the commit message. +env: + has_testthat: 'true' + run_covr: 'true' + run_pkgdown: 'true' + has_RUnit: 'false' + cache-version: 'cache-v1' + run_docker: 'false' + +jobs: + build-check: + runs-on: ${{ matrix.config.os }} + name: ${{ matrix.config.os }} (${{ matrix.config.r }}) + container: ${{ matrix.config.cont }} + ## Environment variables unique to this job. + + strategy: + fail-fast: false + matrix: + config: + - { os: ubuntu-latest, r: 'devel', bioc: '3.17', cont: "bioconductor/bioconductor_docker:RELEASE_3_17", rspm: "https://packagemanager.rstudio.com/cran/__linux__/focal/latest" } + #- { os: macOS-latest, r: 'devel', bioc: '3.15'} + #- { os: windows-latest, r: 'devel', bioc: '3.15'} + # - { os: ubuntu-latest, r: '4.1', bioc: '3.14', cont: "bioconductor/bioconductor_docker:RELEASE_3_14", rspm: "https://packagemanager.rstudio.com/cran/__linux__/focal/latest" } + #- { os: macOS-latest, r: '4.1', bioc: '3.14'} + #- { os: windows-latest, r: '4.1', bioc: '3.14'} + env: + R_REMOTES_NO_ERRORS_FROM_WARNINGS: true + RSPM: ${{ matrix.config.rspm }} + NOT_CRAN: true + TZ: UTC + GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }} + GITHUB_PAT: ${{ secrets.GITHUB_TOKEN }} + + steps: + + ## Set the R library to the directory matching the + ## R packages cache step further below when running on Docker (Linux). + - name: Set R Library home on Linux + if: runner.os == 'Linux' + run: | + mkdir /__w/_temp/Library + echo ".libPaths('/__w/_temp/Library')" > ~/.Rprofile + ## Most of these steps are the same as the ones in + ## https://github.com/r-lib/actions/blob/master/examples/check-standard.yaml + ## If they update their steps, we will also need to update ours. + - name: Checkout Repository + uses: actions/checkout@v2 + + ## R is already included in the Bioconductor docker images + - name: Setup R from r-lib + if: runner.os != 'Linux' + uses: r-lib/actions/setup-r@v2-branch + with: + r-version: ${{ matrix.config.r }} + + ## pandoc is already included in the Bioconductor docker images + - name: Setup pandoc from r-lib + if: runner.os != 'Linux' + uses: r-lib/actions/setup-pandoc@v2-branch + + - name: Query dependencies + run: | + install.packages('remotes') + saveRDS(remotes::dev_package_deps(dependencies = TRUE), ".github/depends.Rds", version = 2) + shell: Rscript {0} + + - name: Cache R packages + if: "!contains(github.event.head_commit.message, '/nocache') && runner.os != 'Linux'" + uses: actions/cache@v2 + with: + path: ${{ env.R_LIBS_USER }} + key: ${{ env.cache-version }}-${{ runner.os }}-biocversion-${{ matrix.config.bioc }}-r-${{ matrix.config.r }}-${{ hashFiles('.github/depends.Rds') }} + restore-keys: ${{ env.cache-version }}-${{ runner.os }}-biocversion-${{ matrix.config.bioc }}-r-${{ matrix.config.r }}-${{ hashFiles('.github/depends.Rds') }} + + - name: Cache R packages on Linux + if: "!contains(github.event.head_commit.message, '/nocache') && runner.os == 'Linux' " + uses: actions/cache@v2 + with: + path: /home/runner/work/_temp/Library + key: ${{ env.cache-version }}-${{ runner.os }}-biocversion-${{ matrix.config.bioc }}-r-${{ matrix.config.r }}-${{ hashFiles('.github/depends.Rds') }} + restore-keys: ${{ env.cache-version }}-${{ runner.os }}-biocversion-${{ matrix.config.bioc }}-r-${{ matrix.config.r }}-${{ hashFiles('.github/depends.Rds') }} + + - name: Install Linux system dependencies + if: runner.os == 'Linux' + run: | + sysreqs=$(Rscript -e 'cat("apt-get update -y && apt-get install -y", paste(gsub("apt-get install -y ", "", remotes::system_requirements("ubuntu", "20.04")), collapse = " "))') + echo $sysreqs + sudo -s eval "$sysreqs" + - name: Install macOS system dependencies + if: matrix.config.os == 'macOS-latest' + run: | + ## Enable installing XML from source if needed + brew install libxml2 + echo "XML_CONFIG=/usr/local/opt/libxml2/bin/xml2-config" >> $GITHUB_ENV + ## Required to install magick as noted at + ## https://github.com/r-lib/usethis/commit/f1f1e0d10c1ebc75fd4c18fa7e2de4551fd9978f#diff-9bfee71065492f63457918efcd912cf2 + brew install imagemagick@6 + ## For textshaping, required by ragg, and required by pkgdown + brew install harfbuzz fribidi + ## For installing usethis's dependency gert + brew install libgit2 + - name: Install Windows system dependencies + if: runner.os == 'Windows' + run: | + ## Edit below if you have any Windows system dependencies + shell: Rscript {0} + + - name: Install BiocManager + run: | + message(paste('****', Sys.time(), 'installing BiocManager ****')) + remotes::install_cran("BiocManager") + shell: Rscript {0} + + - name: Set BiocVersion + run: | + BiocManager::install(version = "${{ matrix.config.bioc }}", ask = FALSE) + shell: Rscript {0} + + - name: Install dependencies pass 1 + run: | + ## Try installing the package dependencies in steps. First the local + ## dependencies, then any remaining dependencies to avoid the + ## issues described at + ## https://stat.ethz.ch/pipermail/bioc-devel/2020-April/016675.html + ## https://github.com/r-lib/remotes/issues/296 + ## Ideally, all dependencies should get installed in the first pass. + ## Pass #1 at installing dependencies + message(paste('****', Sys.time(), 'pass number 1 at installing dependencies: local dependencies ****')) + remotes::install_local(dependencies = TRUE, repos = BiocManager::repositories(), build_vignettes = TRUE, upgrade = TRUE) + continue-on-error: true + shell: Rscript {0} + + - name: Install dependencies pass 2 + run: | + ## Pass #2 at installing dependencies + message(paste('****', Sys.time(), 'pass number 2 at installing dependencies: any remaining dependencies ****')) + remotes::install_local(dependencies = TRUE, repos = BiocManager::repositories(), build_vignettes = TRUE, upgrade = TRUE) + + ## For running the checks + message(paste('****', Sys.time(), 'installing rcmdcheck and BiocCheck ****')) + remotes::install_cran("rcmdcheck") + BiocManager::install("BiocCheck") + shell: Rscript {0} + + - name: Install BiocGenerics + if: env.has_RUnit == 'true' + run: | + ## Install BiocGenerics + BiocManager::install("BiocGenerics") + shell: Rscript {0} + + - name: Install covr + if: github.ref == 'refs/heads/devel' && env.run_covr == 'true' && runner.os == 'Linux' && matrix.config.r == 'devel' + run: | + remotes::install_cran("covr") + shell: Rscript {0} + + - name: Install pkgdown and deps + if: github.ref == 'refs/heads/devel' && env.run_pkgdown == 'true' && runner.os == 'Linux' && matrix.config.r == 'devel' + run: | + remotes::install_cran("pkgdown") + shell: Rscript {0} + + - name: Session info + run: | + options(width = 100) + pkgs <- installed.packages()[, "Package"] + sessioninfo::session_info(pkgs, include_base = TRUE) + shell: Rscript {0} + + - name: Run CMD check + env: + _R_CHECK_CRAN_INCOMING_: false + run: | + rcmdcheck::rcmdcheck( + args = c("--no-build-vignettes", "--no-manual", "--timings"), + build_args = c("--no-manual", "--no-resave-data"), + error_on = "warning", + check_dir = "check" + ) + shell: Rscript {0} + + ## Might need an to add this to the if: && runner.os == 'Linux' + - name: Reveal testthat details + if: env.has_testthat == 'true' + run: find . -name testthat.Rout -exec cat '{}' ';' + + - name: Run RUnit tests + if: env.has_RUnit == 'true' + run: | + BiocGenerics:::testPackage() + shell: Rscript {0} + + - name: Run BiocCheck + env: + DISPLAY: 99.0 + run: | + BiocCheck::BiocCheck( + dir('check', 'tar.gz$', full.names = TRUE), + `quit-with-status` = TRUE, + `no-check-R-ver` = TRUE, + `no-check-bioc-help` = TRUE + ) + shell: Rscript {0} + + - name: Test coverage + if: github.ref == 'refs/heads/devel' && env.run_covr == 'true' && runner.os == 'Linux' && matrix.config.r == 'devel' + run: | + covr::codecov(type="all", commentDontrun=FALSE) + shell: Rscript {0} + + - name: Install package + if: github.ref == 'refs/heads/devel' && env.run_pkgdown == 'true' && runner.os == 'Linux' && matrix.config.r == 'devel' + run: R CMD INSTALL . + + - name: Deploy package + if: github.ref == 'refs/heads/master' && env.run_pkgdown == 'true' && runner.os == 'Linux' && matrix.config.r == 'devel' + run: | + git config --local user.name "$GITHUB_ACTOR" + git config --local user.email "$GITHUB_ACTOR@users.noreply.github.com" + Rscript -e "pkgdown::deploy_to_branch(new_process = FALSE)" + shell: bash {0} + ## Note that you need to run pkgdown::deploy_to_branch(new_process = FALSE) + ## at least one locally before this will work. This creates the gh-pages + ## branch (erasing anything you haven't version controlled!) and + ## makes the git history recognizable by pkgdown. + + - name: Upload check results + if: failure() + uses: actions/upload-artifact@master + with: + name: ${{ runner.os }}-biocversion-${{ matrix.config.bioc }}-r-${{ matrix.config.r }}-results + path: check + + - uses: docker/build-push-action@v1 + if: "!contains(github.event.head_commit.message, '/nodocker') && env.run_docker == 'true' && runner.os == 'Linux' " + with: + username: ${{ secrets.DOCKER_USERNAME }} + password: ${{ secrets.DOCKER_PASSWORD }} + repository: davislaboratory/standR + tag_with_ref: true + tag_with_sha: true + tags: latest diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..816096c --- /dev/null +++ b/.gitignore @@ -0,0 +1,46 @@ +# History files +.Rhistory +.Rapp.history +.Rproj.user +*.Rproj + +# Session Data files +.RData + +# User-specific files +.Ruserdata + +# Example code in package build process +*-Ex.R + +# Output files from R CMD build +/*.tar.gz + +# Output files from R CMD check +/*.Rcheck/ + +# RStudio files +.Rproj.user/ + +# produced vignettes +vignettes/*.html +vignettes/*.pdf + +# OAuth2 token, see https://github.com/hadley/httr/releases/tag/v0.3 +.httr-oauth + +# knitr and R markdown default cache directories +*_cache/ +/cache/ + +# Temporary files created by R markdown +*.utf8.md +*.knit.md + +# R Environment Variables +.Renviron +inst/doc +doc + +# Rproj +hoodscanR.Rproj diff --git a/DESCRIPTION b/DESCRIPTION new file mode 100644 index 0000000..13cfa75 --- /dev/null +++ b/DESCRIPTION @@ -0,0 +1,56 @@ +Package: hoodscanR +Title: Spatial cellular neighbourhood scanning in R +Version: 0.99.0 +Authors@R: + c(person(given = "Ning", + family = "Liu", + role = c("aut", "cre"), + email = "liu.n@wehi.edu.au", + comment = c(ORCID = "0000-0002-9487-9305")), + person(given = "Jarryd", + family = "Martin", + role = c("aut"), + email = "", + comment = c(ORCID = ""))) +Description: hoodscanR is an user-friendly R package providing functions to assist + cellular neighborhood analysis of any spatial transcriptomics data with single-cell resolution. + All functions in the package are built based on the SpatialExperiment object, + allowing integration into various spatial transcriptomics-related packages from Bioconductor. + The package can result in cell-level neighborhood annotation output, along with funtions to + perform neighborhood colocalization analysis and neighborhood-based cell clustering. +biocViews: Spatial, Transcriptomics, SingleCell, Clustering +License: GPL-3 + file LICENSE +URL: https://github.com/DavisLaboratory/hoodscanR +BugReports: https://github.com/DavisLaboratory/hoodscanR/issues +Encoding: UTF-8 +LazyData: false +Roxygen: list(markdown = TRUE) +RoxygenNote: 7.2.3 +Imports: + knitr, + rmarkdown, + SpatialExperiment, + SummarizedExperiment, + circlize, + ComplexHeatmap, + scico, + rlang, + S4Vectors, + dplyr, + tibble, + utils, + ggplot2, + tidyr, + grid, + methods, + stats, + RANN, + Rcpp (>= 1.0.9) +LinkingTo: + Rcpp +Suggests: + testthat (>= 3.0.0) +Config/testthat/edition: 3 +Depends: + R (>= 4.3) +VignetteBuilder: knitr diff --git a/LICENSE b/LICENSE new file mode 100644 index 0000000..94a9ed0 --- /dev/null +++ b/LICENSE @@ -0,0 +1,674 @@ + GNU GENERAL PUBLIC LICENSE + Version 3, 29 June 2007 + + Copyright (C) 2007 Free Software Foundation, Inc. + Everyone is permitted to copy and distribute verbatim copies + of this license document, but changing it is not allowed. + + Preamble + + The GNU General Public License is a free, copyleft license for +software and other kinds of works. + + The licenses for most software and other practical works are designed +to take away your freedom to share and change the works. By contrast, +the GNU General Public License is intended to guarantee your freedom to +share and change all versions of a program--to make sure it remains free +software for all its users. We, the Free Software Foundation, use the +GNU General Public License for most of our software; it applies also to +any other work released this way by its authors. You can apply it to +your programs, too. + + When we speak of free software, we are referring to freedom, not +price. Our General Public Licenses are designed to make sure that you +have the freedom to distribute copies of free software (and charge for +them if you wish), that you receive source code or can get it if you +want it, that you can change the software or use pieces of it in new +free programs, and that you know you can do these things. + + To protect your rights, we need to prevent others from denying you +these rights or asking you to surrender the rights. Therefore, you have +certain responsibilities if you distribute copies of the software, or if +you modify it: responsibilities to respect the freedom of others. + + For example, if you distribute copies of such a program, whether +gratis or for a fee, you must pass on to the recipients the same +freedoms that you received. You must make sure that they, too, receive +or can get the source code. And you must show them these terms so they +know their rights. + + Developers that use the GNU GPL protect your rights with two steps: +(1) assert copyright on the software, and (2) offer you this License +giving you legal permission to copy, distribute and/or modify it. + + For the developers' and authors' protection, the GPL clearly explains +that there is no warranty for this free software. For both users' and +authors' sake, the GPL requires that modified versions be marked as +changed, so that their problems will not be attributed erroneously to +authors of previous versions. + + Some devices are designed to deny users access to install or run +modified versions of the software inside them, although the manufacturer +can do so. This is fundamentally incompatible with the aim of +protecting users' freedom to change the software. The systematic +pattern of such abuse occurs in the area of products for individuals to +use, which is precisely where it is most unacceptable. Therefore, we +have designed this version of the GPL to prohibit the practice for those +products. If such problems arise substantially in other domains, we +stand ready to extend this provision to those domains in future versions +of the GPL, as needed to protect the freedom of users. + + Finally, every program is threatened constantly by software patents. +States should not allow patents to restrict development and use of +software on general-purpose computers, but in those that do, we wish to +avoid the special danger that patents applied to a free program could +make it effectively proprietary. To prevent this, the GPL assures that +patents cannot be used to render the program non-free. + + The precise terms and conditions for copying, distribution and +modification follow. + + TERMS AND CONDITIONS + + 0. Definitions. + + "This License" refers to version 3 of the GNU General Public License. + + "Copyright" also means copyright-like laws that apply to other kinds of +works, such as semiconductor masks. + + "The Program" refers to any copyrightable work licensed under this +License. Each licensee is addressed as "you". "Licensees" and +"recipients" may be individuals or organizations. + + To "modify" a work means to copy from or adapt all or part of the work +in a fashion requiring copyright permission, other than the making of an +exact copy. The resulting work is called a "modified version" of the +earlier work or a work "based on" the earlier work. + + A "covered work" means either the unmodified Program or a work based +on the Program. + + To "propagate" a work means to do anything with it that, without +permission, would make you directly or secondarily liable for +infringement under applicable copyright law, except executing it on a +computer or modifying a private copy. Propagation includes copying, +distribution (with or without modification), making available to the +public, and in some countries other activities as well. + + To "convey" a work means any kind of propagation that enables other +parties to make or receive copies. Mere interaction with a user through +a computer network, with no transfer of a copy, is not conveying. + + An interactive user interface displays "Appropriate Legal Notices" +to the extent that it includes a convenient and prominently visible +feature that (1) displays an appropriate copyright notice, and (2) +tells the user that there is no warranty for the work (except to the +extent that warranties are provided), that licensees may convey the +work under this License, and how to view a copy of this License. If +the interface presents a list of user commands or options, such as a +menu, a prominent item in the list meets this criterion. + + 1. Source Code. + + The "source code" for a work means the preferred form of the work +for making modifications to it. "Object code" means any non-source +form of a work. + + A "Standard Interface" means an interface that either is an official +standard defined by a recognized standards body, or, in the case of +interfaces specified for a particular programming language, one that +is widely used among developers working in that language. + + The "System Libraries" of an executable work include anything, other +than the work as a whole, that (a) is included in the normal form of +packaging a Major Component, but which is not part of that Major +Component, and (b) serves only to enable use of the work with that +Major Component, or to implement a Standard Interface for which an +implementation is available to the public in source code form. A +"Major Component", in this context, means a major essential component +(kernel, window system, and so on) of the specific operating system +(if any) on which the executable work runs, or a compiler used to +produce the work, or an object code interpreter used to run it. + + The "Corresponding Source" for a work in object code form means all +the source code needed to generate, install, and (for an executable +work) run the object code and to modify the work, including scripts to +control those activities. However, it does not include the work's +System Libraries, or general-purpose tools or generally available free +programs which are used unmodified in performing those activities but +which are not part of the work. For example, Corresponding Source +includes interface definition files associated with source files for +the work, and the source code for shared libraries and dynamically +linked subprograms that the work is specifically designed to require, +such as by intimate data communication or control flow between those +subprograms and other parts of the work. + + The Corresponding Source need not include anything that users +can regenerate automatically from other parts of the Corresponding +Source. + + The Corresponding Source for a work in source code form is that +same work. + + 2. Basic Permissions. + + All rights granted under this License are granted for the term of +copyright on the Program, and are irrevocable provided the stated +conditions are met. This License explicitly affirms your unlimited +permission to run the unmodified Program. The output from running a +covered work is covered by this License only if the output, given its +content, constitutes a covered work. This License acknowledges your +rights of fair use or other equivalent, as provided by copyright law. + + You may make, run and propagate covered works that you do not +convey, without conditions so long as your license otherwise remains +in force. You may convey covered works to others for the sole purpose +of having them make modifications exclusively for you, or provide you +with facilities for running those works, provided that you comply with +the terms of this License in conveying all material for which you do +not control copyright. Those thus making or running the covered works +for you must do so exclusively on your behalf, under your direction +and control, on terms that prohibit them from making any copies of +your copyrighted material outside their relationship with you. + + Conveying under any other circumstances is permitted solely under +the conditions stated below. Sublicensing is not allowed; section 10 +makes it unnecessary. + + 3. Protecting Users' Legal Rights From Anti-Circumvention Law. + + No covered work shall be deemed part of an effective technological +measure under any applicable law fulfilling obligations under article +11 of the WIPO copyright treaty adopted on 20 December 1996, or +similar laws prohibiting or restricting circumvention of such +measures. + + When you convey a covered work, you waive any legal power to forbid +circumvention of technological measures to the extent such circumvention +is effected by exercising rights under this License with respect to +the covered work, and you disclaim any intention to limit operation or +modification of the work as a means of enforcing, against the work's +users, your or third parties' legal rights to forbid circumvention of +technological measures. + + 4. Conveying Verbatim Copies. + + You may convey verbatim copies of the Program's source code as you +receive it, in any medium, provided that you conspicuously and +appropriately publish on each copy an appropriate copyright notice; +keep intact all notices stating that this License and any +non-permissive terms added in accord with section 7 apply to the code; +keep intact all notices of the absence of any warranty; and give all +recipients a copy of this License along with the Program. + + You may charge any price or no price for each copy that you convey, +and you may offer support or warranty protection for a fee. + + 5. Conveying Modified Source Versions. + + You may convey a work based on the Program, or the modifications to +produce it from the Program, in the form of source code under the +terms of section 4, provided that you also meet all of these conditions: + + a) The work must carry prominent notices stating that you modified + it, and giving a relevant date. + + b) The work must carry prominent notices stating that it is + released under this License and any conditions added under section + 7. This requirement modifies the requirement in section 4 to + "keep intact all notices". + + c) You must license the entire work, as a whole, under this + License to anyone who comes into possession of a copy. This + License will therefore apply, along with any applicable section 7 + additional terms, to the whole of the work, and all its parts, + regardless of how they are packaged. This License gives no + permission to license the work in any other way, but it does not + invalidate such permission if you have separately received it. + + d) If the work has interactive user interfaces, each must display + Appropriate Legal Notices; however, if the Program has interactive + interfaces that do not display Appropriate Legal Notices, your + work need not make them do so. + + A compilation of a covered work with other separate and independent +works, which are not by their nature extensions of the covered work, +and which are not combined with it such as to form a larger program, +in or on a volume of a storage or distribution medium, is called an +"aggregate" if the compilation and its resulting copyright are not +used to limit the access or legal rights of the compilation's users +beyond what the individual works permit. Inclusion of a covered work +in an aggregate does not cause this License to apply to the other +parts of the aggregate. + + 6. Conveying Non-Source Forms. + + You may convey a covered work in object code form under the terms +of sections 4 and 5, provided that you also convey the +machine-readable Corresponding Source under the terms of this License, +in one of these ways: + + a) Convey the object code in, or embodied in, a physical product + (including a physical distribution medium), accompanied by the + Corresponding Source fixed on a durable physical medium + customarily used for software interchange. + + b) Convey the object code in, or embodied in, a physical product + (including a physical distribution medium), accompanied by a + written offer, valid for at least three years and valid for as + long as you offer spare parts or customer support for that product + model, to give anyone who possesses the object code either (1) a + copy of the Corresponding Source for all the software in the + product that is covered by this License, on a durable physical + medium customarily used for software interchange, for a price no + more than your reasonable cost of physically performing this + conveying of source, or (2) access to copy the + Corresponding Source from a network server at no charge. + + c) Convey individual copies of the object code with a copy of the + written offer to provide the Corresponding Source. This + alternative is allowed only occasionally and noncommercially, and + only if you received the object code with such an offer, in accord + with subsection 6b. + + d) Convey the object code by offering access from a designated + place (gratis or for a charge), and offer equivalent access to the + Corresponding Source in the same way through the same place at no + further charge. You need not require recipients to copy the + Corresponding Source along with the object code. If the place to + copy the object code is a network server, the Corresponding Source + may be on a different server (operated by you or a third party) + that supports equivalent copying facilities, provided you maintain + clear directions next to the object code saying where to find the + Corresponding Source. Regardless of what server hosts the + Corresponding Source, you remain obligated to ensure that it is + available for as long as needed to satisfy these requirements. + + e) Convey the object code using peer-to-peer transmission, provided + you inform other peers where the object code and Corresponding + Source of the work are being offered to the general public at no + charge under subsection 6d. + + A separable portion of the object code, whose source code is excluded +from the Corresponding Source as a System Library, need not be +included in conveying the object code work. + + A "User Product" is either (1) a "consumer product", which means any +tangible personal property which is normally used for personal, family, +or household purposes, or (2) anything designed or sold for incorporation +into a dwelling. In determining whether a product is a consumer product, +doubtful cases shall be resolved in favor of coverage. For a particular +product received by a particular user, "normally used" refers to a +typical or common use of that class of product, regardless of the status +of the particular user or of the way in which the particular user +actually uses, or expects or is expected to use, the product. A product +is a consumer product regardless of whether the product has substantial +commercial, industrial or non-consumer uses, unless such uses represent +the only significant mode of use of the product. + + "Installation Information" for a User Product means any methods, +procedures, authorization keys, or other information required to install +and execute modified versions of a covered work in that User Product from +a modified version of its Corresponding Source. The information must +suffice to ensure that the continued functioning of the modified object +code is in no case prevented or interfered with solely because +modification has been made. + + If you convey an object code work under this section in, or with, or +specifically for use in, a User Product, and the conveying occurs as +part of a transaction in which the right of possession and use of the +User Product is transferred to the recipient in perpetuity or for a +fixed term (regardless of how the transaction is characterized), the +Corresponding Source conveyed under this section must be accompanied +by the Installation Information. But this requirement does not apply +if neither you nor any third party retains the ability to install +modified object code on the User Product (for example, the work has +been installed in ROM). + + The requirement to provide Installation Information does not include a +requirement to continue to provide support service, warranty, or updates +for a work that has been modified or installed by the recipient, or for +the User Product in which it has been modified or installed. Access to a +network may be denied when the modification itself materially and +adversely affects the operation of the network or violates the rules and +protocols for communication across the network. + + Corresponding Source conveyed, and Installation Information provided, +in accord with this section must be in a format that is publicly +documented (and with an implementation available to the public in +source code form), and must require no special password or key for +unpacking, reading or copying. + + 7. Additional Terms. + + "Additional permissions" are terms that supplement the terms of this +License by making exceptions from one or more of its conditions. +Additional permissions that are applicable to the entire Program shall +be treated as though they were included in this License, to the extent +that they are valid under applicable law. If additional permissions +apply only to part of the Program, that part may be used separately +under those permissions, but the entire Program remains governed by +this License without regard to the additional permissions. + + When you convey a copy of a covered work, you may at your option +remove any additional permissions from that copy, or from any part of +it. (Additional permissions may be written to require their own +removal in certain cases when you modify the work.) You may place +additional permissions on material, added by you to a covered work, +for which you have or can give appropriate copyright permission. + + Notwithstanding any other provision of this License, for material you +add to a covered work, you may (if authorized by the copyright holders of +that material) supplement the terms of this License with terms: + + a) Disclaiming warranty or limiting liability differently from the + terms of sections 15 and 16 of this License; or + + b) Requiring preservation of specified reasonable legal notices or + author attributions in that material or in the Appropriate Legal + Notices displayed by works containing it; or + + c) Prohibiting misrepresentation of the origin of that material, or + requiring that modified versions of such material be marked in + reasonable ways as different from the original version; or + + d) Limiting the use for publicity purposes of names of licensors or + authors of the material; or + + e) Declining to grant rights under trademark law for use of some + trade names, trademarks, or service marks; or + + f) Requiring indemnification of licensors and authors of that + material by anyone who conveys the material (or modified versions of + it) with contractual assumptions of liability to the recipient, for + any liability that these contractual assumptions directly impose on + those licensors and authors. + + All other non-permissive additional terms are considered "further +restrictions" within the meaning of section 10. If the Program as you +received it, or any part of it, contains a notice stating that it is +governed by this License along with a term that is a further +restriction, you may remove that term. If a license document contains +a further restriction but permits relicensing or conveying under this +License, you may add to a covered work material governed by the terms +of that license document, provided that the further restriction does +not survive such relicensing or conveying. + + If you add terms to a covered work in accord with this section, you +must place, in the relevant source files, a statement of the +additional terms that apply to those files, or a notice indicating +where to find the applicable terms. + + Additional terms, permissive or non-permissive, may be stated in the +form of a separately written license, or stated as exceptions; +the above requirements apply either way. + + 8. Termination. + + You may not propagate or modify a covered work except as expressly +provided under this License. Any attempt otherwise to propagate or +modify it is void, and will automatically terminate your rights under +this License (including any patent licenses granted under the third +paragraph of section 11). + + However, if you cease all violation of this License, then your +license from a particular copyright holder is reinstated (a) +provisionally, unless and until the copyright holder explicitly and +finally terminates your license, and (b) permanently, if the copyright +holder fails to notify you of the violation by some reasonable means +prior to 60 days after the cessation. + + Moreover, your license from a particular copyright holder is +reinstated permanently if the copyright holder notifies you of the +violation by some reasonable means, this is the first time you have +received notice of violation of this License (for any work) from that +copyright holder, and you cure the violation prior to 30 days after +your receipt of the notice. + + Termination of your rights under this section does not terminate the +licenses of parties who have received copies or rights from you under +this License. If your rights have been terminated and not permanently +reinstated, you do not qualify to receive new licenses for the same +material under section 10. + + 9. Acceptance Not Required for Having Copies. + + You are not required to accept this License in order to receive or +run a copy of the Program. Ancillary propagation of a covered work +occurring solely as a consequence of using peer-to-peer transmission +to receive a copy likewise does not require acceptance. However, +nothing other than this License grants you permission to propagate or +modify any covered work. These actions infringe copyright if you do +not accept this License. Therefore, by modifying or propagating a +covered work, you indicate your acceptance of this License to do so. + + 10. Automatic Licensing of Downstream Recipients. + + Each time you convey a covered work, the recipient automatically +receives a license from the original licensors, to run, modify and +propagate that work, subject to this License. You are not responsible +for enforcing compliance by third parties with this License. + + An "entity transaction" is a transaction transferring control of an +organization, or substantially all assets of one, or subdividing an +organization, or merging organizations. If propagation of a covered +work results from an entity transaction, each party to that +transaction who receives a copy of the work also receives whatever +licenses to the work the party's predecessor in interest had or could +give under the previous paragraph, plus a right to possession of the +Corresponding Source of the work from the predecessor in interest, if +the predecessor has it or can get it with reasonable efforts. + + You may not impose any further restrictions on the exercise of the +rights granted or affirmed under this License. For example, you may +not impose a license fee, royalty, or other charge for exercise of +rights granted under this License, and you may not initiate litigation +(including a cross-claim or counterclaim in a lawsuit) alleging that +any patent claim is infringed by making, using, selling, offering for +sale, or importing the Program or any portion of it. + + 11. Patents. + + A "contributor" is a copyright holder who authorizes use under this +License of the Program or a work on which the Program is based. The +work thus licensed is called the contributor's "contributor version". + + A contributor's "essential patent claims" are all patent claims +owned or controlled by the contributor, whether already acquired or +hereafter acquired, that would be infringed by some manner, permitted +by this License, of making, using, or selling its contributor version, +but do not include claims that would be infringed only as a +consequence of further modification of the contributor version. For +purposes of this definition, "control" includes the right to grant +patent sublicenses in a manner consistent with the requirements of +this License. + + Each contributor grants you a non-exclusive, worldwide, royalty-free +patent license under the contributor's essential patent claims, to +make, use, sell, offer for sale, import and otherwise run, modify and +propagate the contents of its contributor version. + + In the following three paragraphs, a "patent license" is any express +agreement or commitment, however denominated, not to enforce a patent +(such as an express permission to practice a patent or covenant not to +sue for patent infringement). To "grant" such a patent license to a +party means to make such an agreement or commitment not to enforce a +patent against the party. + + If you convey a covered work, knowingly relying on a patent license, +and the Corresponding Source of the work is not available for anyone +to copy, free of charge and under the terms of this License, through a +publicly available network server or other readily accessible means, +then you must either (1) cause the Corresponding Source to be so +available, or (2) arrange to deprive yourself of the benefit of the +patent license for this particular work, or (3) arrange, in a manner +consistent with the requirements of this License, to extend the patent +license to downstream recipients. "Knowingly relying" means you have +actual knowledge that, but for the patent license, your conveying the +covered work in a country, or your recipient's use of the covered work +in a country, would infringe one or more identifiable patents in that +country that you have reason to believe are valid. + + If, pursuant to or in connection with a single transaction or +arrangement, you convey, or propagate by procuring conveyance of, a +covered work, and grant a patent license to some of the parties +receiving the covered work authorizing them to use, propagate, modify +or convey a specific copy of the covered work, then the patent license +you grant is automatically extended to all recipients of the covered +work and works based on it. + + A patent license is "discriminatory" if it does not include within +the scope of its coverage, prohibits the exercise of, or is +conditioned on the non-exercise of one or more of the rights that are +specifically granted under this License. You may not convey a covered +work if you are a party to an arrangement with a third party that is +in the business of distributing software, under which you make payment +to the third party based on the extent of your activity of conveying +the work, and under which the third party grants, to any of the +parties who would receive the covered work from you, a discriminatory +patent license (a) in connection with copies of the covered work +conveyed by you (or copies made from those copies), or (b) primarily +for and in connection with specific products or compilations that +contain the covered work, unless you entered into that arrangement, +or that patent license was granted, prior to 28 March 2007. + + Nothing in this License shall be construed as excluding or limiting +any implied license or other defenses to infringement that may +otherwise be available to you under applicable patent law. + + 12. No Surrender of Others' Freedom. + + If conditions are imposed on you (whether by court order, agreement or +otherwise) that contradict the conditions of this License, they do not +excuse you from the conditions of this License. If you cannot convey a +covered work so as to satisfy simultaneously your obligations under this +License and any other pertinent obligations, then as a consequence you may +not convey it at all. For example, if you agree to terms that obligate you +to collect a royalty for further conveying from those to whom you convey +the Program, the only way you could satisfy both those terms and this +License would be to refrain entirely from conveying the Program. + + 13. Use with the GNU Affero General Public License. + + Notwithstanding any other provision of this License, you have +permission to link or combine any covered work with a work licensed +under version 3 of the GNU Affero General Public License into a single +combined work, and to convey the resulting work. The terms of this +License will continue to apply to the part which is the covered work, +but the special requirements of the GNU Affero General Public License, +section 13, concerning interaction through a network will apply to the +combination as such. + + 14. Revised Versions of this License. + + The Free Software Foundation may publish revised and/or new versions of +the GNU General Public License from time to time. Such new versions will +be similar in spirit to the present version, but may differ in detail to +address new problems or concerns. + + Each version is given a distinguishing version number. If the +Program specifies that a certain numbered version of the GNU General +Public License "or any later version" applies to it, you have the +option of following the terms and conditions either of that numbered +version or of any later version published by the Free Software +Foundation. If the Program does not specify a version number of the +GNU General Public License, you may choose any version ever published +by the Free Software Foundation. + + If the Program specifies that a proxy can decide which future +versions of the GNU General Public License can be used, that proxy's +public statement of acceptance of a version permanently authorizes you +to choose that version for the Program. + + Later license versions may give you additional or different +permissions. However, no additional obligations are imposed on any +author or copyright holder as a result of your choosing to follow a +later version. + + 15. Disclaimer of Warranty. + + THERE IS NO WARRANTY FOR THE PROGRAM, TO THE EXTENT PERMITTED BY +APPLICABLE LAW. EXCEPT WHEN OTHERWISE STATED IN WRITING THE COPYRIGHT +HOLDERS AND/OR OTHER PARTIES PROVIDE THE PROGRAM "AS IS" WITHOUT WARRANTY +OF ANY KIND, EITHER EXPRESSED OR IMPLIED, INCLUDING, BUT NOT LIMITED TO, +THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR +PURPOSE. THE ENTIRE RISK AS TO THE QUALITY AND PERFORMANCE OF THE PROGRAM +IS WITH YOU. SHOULD THE PROGRAM PROVE DEFECTIVE, YOU ASSUME THE COST OF +ALL NECESSARY SERVICING, REPAIR OR CORRECTION. + + 16. Limitation of Liability. + + IN NO EVENT UNLESS REQUIRED BY APPLICABLE LAW OR AGREED TO IN WRITING +WILL ANY COPYRIGHT HOLDER, OR ANY OTHER PARTY WHO MODIFIES AND/OR CONVEYS +THE PROGRAM AS PERMITTED ABOVE, BE LIABLE TO YOU FOR DAMAGES, INCLUDING ANY +GENERAL, SPECIAL, INCIDENTAL OR CONSEQUENTIAL DAMAGES ARISING OUT OF THE +USE OR INABILITY TO USE THE PROGRAM (INCLUDING BUT NOT LIMITED TO LOSS OF +DATA OR DATA BEING RENDERED INACCURATE OR LOSSES SUSTAINED BY YOU OR THIRD +PARTIES OR A FAILURE OF THE PROGRAM TO OPERATE WITH ANY OTHER PROGRAMS), +EVEN IF SUCH HOLDER OR OTHER PARTY HAS BEEN ADVISED OF THE POSSIBILITY OF +SUCH DAMAGES. + + 17. Interpretation of Sections 15 and 16. + + If the disclaimer of warranty and limitation of liability provided +above cannot be given local legal effect according to their terms, +reviewing courts shall apply local law that most closely approximates +an absolute waiver of all civil liability in connection with the +Program, unless a warranty or assumption of liability accompanies a +copy of the Program in return for a fee. + + END OF TERMS AND CONDITIONS + + How to Apply These Terms to Your New Programs + + If you develop a new program, and you want it to be of the greatest +possible use to the public, the best way to achieve this is to make it +free software which everyone can redistribute and change under these terms. + + To do so, attach the following notices to the program. It is safest +to attach them to the start of each source file to most effectively +state the exclusion of warranty; and each file should have at least +the "copyright" line and a pointer to where the full notice is found. + + + Copyright (C) + + This program is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see . + +Also add information on how to contact you by electronic and paper mail. + + If the program does terminal interaction, make it output a short +notice like this when it starts in an interactive mode: + + Copyright (C) + This program comes with ABSOLUTELY NO WARRANTY; for details type `show w'. + This is free software, and you are welcome to redistribute it + under certain conditions; type `show c' for details. + +The hypothetical commands `show w' and `show c' should show the appropriate +parts of the General Public License. Of course, your program's commands +might be different; for a GUI interface, you would use an "about box". + + You should also get your employer (if you work as a programmer) or school, +if any, to sign a "copyright disclaimer" for the program, if necessary. +For more information on this, and how to apply and follow the GNU GPL, see +. + + The GNU General Public License does not permit incorporating your program +into proprietary programs. If your program is a subroutine library, you +may consider it more useful to permit linking proprietary applications with +the library. If this is what you want to do, use the GNU Lesser General +Public License instead of this License. But first, please read +. diff --git a/NAMESPACE b/NAMESPACE new file mode 100644 index 0000000..5907c4c --- /dev/null +++ b/NAMESPACE @@ -0,0 +1,37 @@ +# Generated by roxygen2: do not edit by hand + +export(calcMetrics) +export(clustByHood) +export(findNearCells) +export(mergeByGroup) +export(mergeHoodSpe) +export(plotColocal) +export(plotHoodMat) +export(plotProbDist) +export(plotTissue) +export(readHoodData) +export(scanHoods) +import(Rcpp) +import(SpatialExperiment) +import(ggplot2) +import(knitr) +importFrom(S4Vectors,"metadata<-") +importFrom(S4Vectors,metadata) +importFrom(SummarizedExperiment,assay) +importFrom(SummarizedExperiment,colData) +importFrom(dplyr,across) +importFrom(dplyr,filter) +importFrom(dplyr,group_by) +importFrom(dplyr,summarise) +importFrom(dplyr,where) +importFrom(grid,gpar) +importFrom(grid,unit) +importFrom(methods,is) +importFrom(rmarkdown,draft) +importFrom(stats,cor) +importFrom(stats,median) +importFrom(tibble,column_to_rownames) +importFrom(tibble,rownames_to_column) +importFrom(tidyr,pivot_longer) +importFrom(utils,globalVariables) +useDynLib(hoodscanR) diff --git a/NEWS.md b/NEWS.md new file mode 100644 index 0000000..7d0f2fb --- /dev/null +++ b/NEWS.md @@ -0,0 +1,2 @@ +# hoodscanR 1.0.0 +First release. diff --git a/R/RcppExports.R b/R/RcppExports.R new file mode 100644 index 0000000..51aa637 --- /dev/null +++ b/R/RcppExports.R @@ -0,0 +1,7 @@ +# Generated by using Rcpp::compileAttributes() -> do not edit by hand +# Generator token: 10BE3573-1514-4C36-9D1C-5A225CD40393 + +calculate_metrics <- function(probs) { + .Call('_hoodscanR_calculate_metrics', PACKAGE = 'hoodscanR', probs) +} + diff --git a/R/calc_metrics.R b/R/calc_metrics.R new file mode 100644 index 0000000..daba4d9 --- /dev/null +++ b/R/calc_metrics.R @@ -0,0 +1,47 @@ +#' Calculate metrics for probability matrix +#' +#' @param spe A SpatialExperiment object. +#' @param pm Optional. The probability matrix. +#' @param pm_cols The colnames of probability matrix. This is requires for +#' SpatialExperiment input. Assuming that the probability is stored in the colData. +#' +#' @return A SpatialExperiment object. Calculated entropy and perplexity are +#' saved as columns in the colData of the SpatialExperiment object. +#' @export +#' +#' @examples +#' +#' data("spe_test") +#' +#' spe <- readHoodData(spe, anno_col = "celltypes") +#' +#' fnc <- findNearCells(spe, k = 100) +#' +#' pm <- scanHoods(fnc$distance) +#' +#' pm2 <- mergeByGroup(pm, fnc$cells) +#' +#' spe <- mergeHoodSpe(spe, pm2) +#' +#' spe <- calcMetrics(spe, pm_cols = colnames(pm2)) +#' +calcMetrics <- function(spe, pm = NA, pm_cols = NA) { + if (is(pm, "logical")) { + if (is(pm_cols, "logical")) { + stop("Need to input either the pm or pm_cols parameters.") + } else { + pm <- as.data.frame(colData(spe))[, pm_cols] |> + as.matrix() + } + } else { + pm <- pm + } + + result <- calculate_metrics(pm) + + for (i in colnames(result)) { + colData(spe)[, i] <- result[, i] + } + + return(spe) +} diff --git a/R/find_near_cells.R b/R/find_near_cells.R new file mode 100644 index 0000000..9584378 --- /dev/null +++ b/R/find_near_cells.R @@ -0,0 +1,91 @@ +#' Find the k-th nearest cells for each cell +#' +#' @param dat A SpatialExperiment object, can be generated using function `readHoodData`. +#' @param k The maximum number of nearest cells to compute. +#' @param targetCell Specify the cells to be the target cell for finding nearest cells. +#' @param reportCellID Logical. Set to TRUE to report cell id instead of cell types. +#' @param reportDist Logical. Set to TRUE to report the distance matrix. +#' @param anno_col Character vector. The name of annotation column to use. +#' +#' @return A list includes a data.frame and a matrix, describing the cell types +#' and distances of the k-th nearest cells of each cell. +#' @details The `findNearCells` function uses the `nn2` function from the `RANN` package, +#' which uses the Approximate Near Neighbor (ANN) C++ library. For more infromation on the +#' ANN library please see http://www.cs.umd.edu/~mount/ANN/. +#' +#' @export +#' +#' @examples +#' +#' data("spe_test") +#' +#' spe <- readHoodData(spe, anno_col = "celltypes") +#' +#' fnc <- findNearCells(spe, k = 100) +#' +findNearCells <- function(dat, k = 100, targetCell = FALSE, + reportCellID = FALSE, + reportDist = TRUE, anno_col = 0) { + if (!is(dat, "SpatialExperiment")) { + stop("Please use SpatialExperiment as input.") + } + + cell_pos_dat <- SpatialExperiment::spatialCoords(dat) + colnames(cell_pos_dat) <- c("x", "y") + + if (anno_col == 0) { + anno_col <- "cell_annotation" + } + + cell_annotation <- SummarizedExperiment::colData(dat)[, anno_col] + + query <- targetCell + + if (isFALSE(targetCell)) { + query <- cell_pos_dat + } else { + if (length(targetCell) == 1) { + query <- t(as.matrix(cell_pos_dat[targetCell, ])) + rownames(query) <- targetCell + } else { + query <- cell_pos_dat[targetCell, ] + } + } + + searchcells <- RANN::nn2(data = cell_pos_dat, query = query, + k = k + 1, searchtype = "priority") + + closest <- searchcells[[1]] + + idxcell <- closest[, 1] + idxclosecell <- closest[, seq(from = 2, to = ncol(closest))] + + if (isFALSE(reportCellID)) { + idxcell[] <- rownames(cell_pos_dat) + idxclosecell[] <- cell_annotation[c(idxclosecell)] + } else { + idxcell[] <- rownames(cell_pos_dat)[c(idxcell)] + idxclosecell[] <- rownames(cell_pos_dat)[c(idxclosecell)] + } + + if (!isFALSE(targetCell) & length(targetCell) == 1) { + result <- as.data.frame(t(as.data.frame(idxclosecell))) + rownames(result) <- targetCell + } else { + result <- cbind(idxcell, idxclosecell) |> + as.data.frame() |> + tibble::column_to_rownames("idxcell") + } + + colnames(result) <- paste0("nearest_cell_", seq(k)) + + if (isTRUE(reportDist)) { + dist <- searchcells$nn.dists[, seq(from = 2, to = (k + 1))] + rownames(dist) <- rownames(result) + colnames(dist) <- colnames(result) + result <- list(result, dist) + names(result) <- c("cells", "distance") + } + + return(result) +} diff --git a/R/hoodscanR-package.R b/R/hoodscanR-package.R new file mode 100644 index 0000000..7bbbe22 --- /dev/null +++ b/R/hoodscanR-package.R @@ -0,0 +1,61 @@ +#' @importFrom SummarizedExperiment colData +#' @import SpatialExperiment +#' @importFrom dplyr across +#' @import ggplot2 +#' @importFrom SummarizedExperiment assay +#' @importFrom tibble column_to_rownames +#' @importFrom tibble rownames_to_column +#' @importFrom dplyr group_by +#' @importFrom S4Vectors metadata +#' @importFrom S4Vectors metadata<- +#' @importFrom tidyr pivot_longer +#' @importFrom dplyr filter +#' @importFrom dplyr summarise +#' @importFrom dplyr where +#' @importFrom methods is +#' @importFrom grid gpar +#' @importFrom grid unit +#' @importFrom stats cor +#' @importFrom stats median +#' @import Rcpp +#' @import knitr +#' @importFrom utils globalVariables +#' @importFrom rmarkdown draft +NULL + +#' Method to identify cellular spatial neighbourhood from single cell spatial +#' transcriptomics data. +#' +#' `hoodscanR` implements a novel method to scan for cell neighbourhood from +#' spatial transcriptomics data at single cell level, such as CosMx and MERFISH etc. +#' `hoodscanR` takes the cellular position and cell type annotations +#' as inputs, allowing cellular spatial neighbourhood analysis. +#' +#' +#' +#' @author Ning Liu \email{liu.n@@wehi.edu.au} +#' @name hoodscanR-package +#' @docType package +#' @aliases hoodscanR hoodscanR-package +#' @keywords internal +#' +#' @useDynLib hoodscanR +#' +NULL + + +#' Example test data +#' +#' hoodscanR-package has 1 datasets: \itemize{ +#' \item spe_test An example data +#' } +#' +#' @docType data +#' @name spe_test +#' @usage data("spe_test") +#' @return A SpatialExperiment object +#' @keywords internal +#' @format A SpatialExperiment object +#' @examples +#' data(spe_test) +"spe" diff --git a/R/merge_by_group.R b/R/merge_by_group.R new file mode 100644 index 0000000..92c73f5 --- /dev/null +++ b/R/merge_by_group.R @@ -0,0 +1,47 @@ +#' Merge probability matrix based on annotations +#' +#' @param pm A numeric matrix. Probability matrix generated by the soft_max function. +#' @param group_df A character matrix. Annotation of the neighboring cells to be used. +#' +#' @return A probability matrix, describing the probability of each cell being +#' in each cellular neighborhood. +#' @export +#' +#' @examples +#' +#' data("spe_test") +#' +#' spe <- readHoodData(spe, anno_col = "celltypes") +#' +#' fnc <- findNearCells(spe, k = 100) +#' +#' pm <- scanHoods(fnc$distance) +#' +#' pm2 <- mergeByGroup(pm, fnc$cells) +#' +mergeByGroup <- function(pm, group_df) { + df <- as.matrix(pm) + group_df <- as.matrix(group_df) + + if (!all(dim(df) == dim(group_df))) { + stop("df and group_df should have the same dimensions.") + } + + colnames(df) <- colnames(group_df) + rownames(df) <- rownames(group_df) + + uniquegroups <- sort(unique(as.vector((group_df)))) + + out_df <- matrix(nrow = nrow(pm), ncol = length(uniquegroups)) + colnames(out_df) <- uniquegroups + rownames(out_df) <- rownames(df) + for (g in uniquegroups) { + mask <- (group_df == g) |> `storage.mode<-`("numeric") + + out_df[, g] <- rowSums(mask * pm)[rownames(df)] + } + + out_df[is.na(out_df)] <- 0 + + return(out_df) +} diff --git a/R/merge_pm_spe.R b/R/merge_pm_spe.R new file mode 100644 index 0000000..8585ab7 --- /dev/null +++ b/R/merge_pm_spe.R @@ -0,0 +1,34 @@ +#' Merge probability matrix into SpatialExperiment object. +#' +#' @param spe A SpatialExperiment object. +#' @param pm Probability matrix. Can be obtained by the function mergeByGroup. +#' +#' @return A SpatialExperiment object. Cell-level neighborhood information are +#' saved in the colData of the SpatialExperiment object. +#' @export +#' +#' @examples +#' +#' data("spe_test") +#' +#' spe <- readHoodData(spe, anno_col = "celltypes") +#' +#' fnc <- findNearCells(spe, k = 100) +#' +#' pm <- scanHoods(fnc$distance) +#' +#' pm2 <- mergeByGroup(pm, fnc$cells) +#' +#' spe <- mergeHoodSpe(spe, pm2) +#' +mergeHoodSpe <- function(spe, pm) { + if (!(all(rownames(pm) == colnames(spe)))) { + pm <- pm[colnames(spe), ] + } + + for (i in colnames(pm)) { + colData(spe)[, i] <- pm[, i] + } + + return(spe) +} diff --git a/R/plot_heatmap.R b/R/plot_heatmap.R new file mode 100644 index 0000000..c739494 --- /dev/null +++ b/R/plot_heatmap.R @@ -0,0 +1,150 @@ +#' Plot heatmap for neighbourhood analysis +#' +#' @param object A probability matrix or SpatialExperiment. +#' @param pm_cols The colnames of probability matrix. This is requires for +#' SpatialExperiment input. Assuming that the probability is stored in the colData. +#' @param self.cor Logical. By default is TRUE, inidicating running a correlation +#' between neighbourhoods to perform a simple co-localization analysis. +#' When this set to FALSE, it will plot the average probability of each +#' neighbourhood by group using the byGroup parameter. +#' @param byGroup Character. This is required when self.cor is set to FALSE. +#' @param hm_width Integer. The width of heatmap. +#' @param hm_height Integer. The height of heatmap. +#' @param clusterRows Logical. Cluster rows. +#' @param clusterCols Logical. Cluster columns. +#' @param return_matrix Logical. Export a numeric matrix . +#' @param ... Ignore parameter. +#' +#' @return A ComplexHeatmap plot. When return_matrix is set to TRUE, +#' return a matrix Object. +#' @export +#' +#' +#' @docType methods +#' @rdname plotColocal +#' +#' @aliases plotColocal plotColocal, matrix-method +#' @aliases plotColocal, SpatialExperiment-method plotColocal +#' +#' @examples +#' +#' data("spe_test") +#' +#' spe <- readHoodData(spe, anno_col = "celltypes") +#' +#' fnc <- findNearCells(spe, k = 100) +#' +#' pm <- scanHoods(fnc$distance) +#' +#' pm2 <- mergeByGroup(pm, fnc$cells) +#' +#' spe <- mergeHoodSpe(spe, pm2) +#' +#' plotColocal(spe, pm_cols = colnames(pm2)) +#' +setGeneric( + "plotColocal", + function(object, ...) standardGeneric("plotColocal") +) + + +#' @rdname plotColocal +setMethod( + "plotColocal", + signature("matrix"), + function(object, hm_width = 5, hm_height = 5) { + cor.m <- cor(object) + ht <- plotColocal_intl(cor.m, hm_width = hm_width, hm_height = hm_height) + + ComplexHeatmap::draw(ht) + } +) + +#' @rdname plotColocal +setMethod( + "plotColocal", + signature("SpatialExperiment"), + function(object, pm_cols, self.cor = TRUE, byGroup = NA, hm_width = 5, + hm_height = 5, clusterRows = TRUE, clusterCols = TRUE, + return_matrix = FALSE) { + dat <- as.data.frame(colData(object)) + + if (!all(pm_cols %in% colnames(dat))) { + stop("The pm_cols are not included in the SpatialExperiment.") + } + + if (self.cor) { + cor.m <- cor(dat[, pm_cols]) + ht <- plotColocal_intl(cor.m, + hm_width = hm_width, hm_height = hm_height, + clusterRows = clusterRows, clusterCols = clusterCols + ) + } else { + if (is(byGroup, "logical") | length(byGroup) != 1 | !(byGroup %in% colnames(dat))) { + stop("byGroup should be the columns in colData") + } else { + datx <- dat[, c(pm_cols, byGroup)] + cor.m <- mean_by_group(datx, byGroup) |> + as.data.frame() |> + column_to_rownames(byGroup) |> + as.matrix() + ht <- plotColocal_intl(cor.m, + title = "Mean probability within groups", legend.name = "Prob.", + hm_width = hm_width, hm_height = hm_height, self.cor = FALSE + ) + } + } + + if (isTRUE(return_matrix)) { + return(cor.m) + } else { + ComplexHeatmap::draw(ht) + } + } +) + +mean_by_group <- function(df, group_col) { + df[, group_col] <- as.character(df[, group_col]) + + df_out <- df |> + group_by(!!rlang::sym(group_col)) |> + summarise(across(where(is.numeric), mean)) + + return(df_out) +} + +plotColocal_intl <- function(m, col.pal = NA, + title = "Pearson correlation between neighbourhoods", + title.size = 15, legend.name = "Cor.", + hm_width, hm_height, self.cor = TRUE, + clusterRows = FALSE, clusterCols = FALSE) { + if (isTRUE(self.cor)) { + if (is(col.pal, "logical")) { + cp <- circlize::colorRamp2(seq(-1, 1, by = 0.01), + scico::scico(201, palette = "roma", direction = -1)) + } else { + cp <- col.pal + } + } else { + if (is(col.pal, "logical")) { + cp <- circlize::colorRamp2(seq(0, 1, by = 0.01), + scico::scico(101, palette = "lajolla")) + } else { + cp <- col.pal + } + } + + + + ht <- ComplexHeatmap::Heatmap(m, + col = cp, rect_gp = gpar(col = "white", lwd = 2), + name = legend.name, column_title = title, + column_title_gp = gpar(fontsize = title.size, fontface = "italic"), + column_dend_height = unit(1, "cm"), + row_dend_width = unit(1, "cm"), column_names_rot = 45, + width = unit(hm_width, "cm"), height = unit(hm_height, "cm"), + cluster_rows = clusterRows, cluster_columns = clusterCols + ) + + return(ht) +} diff --git a/R/plot_pd.R b/R/plot_pd.R new file mode 100644 index 0000000..bed0d66 --- /dev/null +++ b/R/plot_pd.R @@ -0,0 +1,199 @@ +#' Plot probability distribution +#' +#' @param object A probability matrix or SpatialExperiment. +#' @param targetCells Character. Optional. Can speicify one or more cells to be plotted. +#' @param pm_cols The colnames of probability matrix. This is requires for +#' SpatialExperiment input. Assuming that the probability is stored in the colData. +#' @param by_cluster Logical. By default is TRUE, to plot distribution by each cluster. +#' @param show_clusters Character. The cluster to be ploted, by default is 1 to 6. +#' @param plot_all Logical. By default is FALSE, set this to true to plot box +#' plot instead of bar plot to show all cells in each cluster. +#' @param sample_size Integer. By default is 2, sampling two cell from each +#' cluster to be plotted. +#' @param ... aesthetic mappings to pass to `ggplot2::aes_string()`. +#' +#' @return A ggplot object. +#' @export +#' +#' @docType methods +#' @name plotProbDist +#' @rdname plotProbDist +#' @aliases plotProbDist plotProbDist,matrix-method plotProbDist,SpatialExperiment-method +#' +#' @examples +#' +#' data("spe_test") +#' +#' spe <- readHoodData(spe, anno_col = "celltypes") +#' +#' fnc <- findNearCells(spe, k = 100) +#' +#' pm <- scanHoods(fnc$distance) +#' +#' pm2 <- mergeByGroup(pm, fnc$cells) +#' +#' spe <- mergeHoodSpe(spe, pm2) +#' +#' plotProbDist(spe, pm_cols = colnames(pm2)) +#' +setGeneric( + "plotProbDist", + function(object, ...) standardGeneric("plotProbDist") +) + + +#' @rdname plotProbDist +setMethod( + "plotProbDist", + signature("matrix"), + function(object, targetCells = NA, ...) { + if (is(targetCells, "logical")) { + targetCells <- rownames(object)[seq(6)] + } + + dat <- as.data.frame(object) |> + rownames_to_column("cells") |> + filter(cells %in% targetCells) |> + pivot_longer(cols = -cells, names_to = "hoods", values_to = "probability") + + p <- plotProbDist_intl(dat, ...) + + facet_wrap(~cells, ncol = 4) + + return(p) + } +) + +#' @rdname plotProbDist +setMethod( + "plotProbDist", + signature("SpatialExperiment"), + function(object, pm_cols, targetCells = NA, by_cluster = FALSE, + show_clusters = as.character(seq(6)), plot_all = FALSE, + sample_size = 2, ...) { + dat <- as.data.frame(colData(object)) + + if (!all(pm_cols %in% colnames(dat))) { + stop("The pm_cols are not included in the SpatialExperiment.") + } + + if (isTRUE(by_cluster)) { + if (!("clusters" %in% colnames(dat))) { + stop("Cannot find the clusters column in the SpatialExperiment, check column names.") + } + + if (all(!(show_clusters %in% dat[, "clusters"]))) { + stop("Cannot find the show_clusters value in the clusters column.") + } + + dat <- dat[, c(pm_cols, "clusters")] + + dat <- as.data.frame(dat) |> + rownames_to_column("cells") |> + filter(clusters %in% show_clusters) + + if (isTRUE(plot_all)) { + p <- plotProbDist_box_intl(dat, pm_cols) + + facet_wrap(~clusters) + } else { + datx <- sample_rows(dat, "clusters", sample_size) + + datx <- datx |> + pivot_longer(cols = -c(cells, clusters), names_to = "hoods", + values_to = "probability") + + p <- plotProbDist_intl(datx, ...) + + facet_wrap(~ clusters + cells, ncol = 4) + } + } else { + dat <- dat[, c(pm_cols)] + + if (is(targetCells, "logical")) { + targetCells <- rownames(dat)[seq(6)] + } + + dat <- as.data.frame(dat) |> + rownames_to_column("cells") |> + filter(cells %in% targetCells) |> + pivot_longer(cols = -cells, names_to = "hoods", + values_to = "probability") + + p <- plotProbDist_intl(dat, ...) + + facet_wrap(~cells, ncol = 4) + } + + return(p) + } +) + + + +sample_rows <- function(df, group_var, sample_size) { + df |> + group_by(!!sym(group_var)) |> + dplyr::slice_sample(n = sample_size, replace = FALSE) +} + + +plotProbDist_intl <- function(x, ...) { + aesmap <- rlang::enquos(...) + + # split aes params into those that are not aes i.e. static parametrisation + if (length(aesmap) > 0) { + is_aes <- vapply(aesmap, rlang::quo_is_symbolic, FUN.VALUE = logical(1)) + defaultmap <- lapply(aesmap[!is_aes], rlang::eval_tidy) + aesmap <- aesmap[is_aes] + } else { + defaultmap <- list() + } + + defaultmap[["stat"]] <- "identity" + + p <- ggplot2::ggplot(x, aes(x = hoods, y = probability, !!!aesmap)) + + p <- p + + do.call(ggplot2::geom_bar, defaultmap) + + pm_theme() + + return(p) +} + +plotProbDist_box_intl <- function(x, pm_cols, ...) { + aesmap <- rlang::enquos(...) + + # split aes params into those that are not aes i.e. static parametrisation + if (length(aesmap) > 0) { + is_aes <- vapply(aesmap, rlang::quo_is_symbolic, FUN.VALUE = logical(1)) + defaultmap <- lapply(aesmap[!is_aes], rlang::eval_tidy) + aesmap <- aesmap[is_aes] + } else { + defaultmap <- list() + } + + x <- pivot_longer(x, + cols = pm_cols, + names_to = "hoods", values_to = "probability" + ) + + p <- ggplot2::ggplot(x, aes(x = hoods, y = probability, !!!aesmap)) + + p <- p + + do.call(ggplot2::geom_boxplot, defaultmap) + + pm_theme() + + return(p) +} + + +pm_theme <- function(textScale = 1.1) { + stopifnot(textScale > 0) + ggplot2::theme_bw() + + ggplot2::theme( + panel.border = element_rect(colour = "black", fill = NA, linewidth = 1), + panel.grid = element_blank(), + legend.text = element_text(size = rel(textScale)), + legend.title = element_text(size = rel(textScale), face = "italic"), + axis.text.x = element_text(angle = 90, vjust = 0.5, hjust = 1) + ) +} + +utils::globalVariables(c("cells", "clusters", "hoods", "probability")) diff --git a/R/plot_pm.R b/R/plot_pm.R new file mode 100644 index 0000000..b88317c --- /dev/null +++ b/R/plot_pm.R @@ -0,0 +1,99 @@ +#' Plot probability matrix as a heatmap +#' +#' @param object A probability matrix or SpatialExperiment. +#' @param pm_cols The colnames of probability matrix. This is requires for +#' SpatialExperiment input. Assuming that the probability is stored in the colData. +#' @param targetCells Character. Optional. Can speicify one or more cells to be plotted. +#' @param n Integer. The number of randomly selected cells to be plotted. +#' This parameter will be used when targetCells is not specify. +#' @param hm_width Integer. The width of heatmap. +#' @param hm_height Integer. The height of heatmap. +#' @param clusterRows Logical. Cluster rows or not. +#' @param clusterCols Logical. Cluster columns or not. +#' @param title Title of the heatmap. +#' @param ... Ignore parameter. +#' +#' @return A ComplexHeatmap plot. +#' @export +#' +#' @docType methods +#' @name plotHoodMat +#' @rdname plotHoodMat +#' @aliases plotHoodMat plotHoodMat,matrix-method plotHoodMat,SpatialExperiment-method +#' +#' @examples +#' +#' data("spe_test") +#' +#' spe <- readHoodData(spe, anno_col = "celltypes") +#' +#' fnc <- findNearCells(spe, k = 100) +#' +#' pm <- scanHoods(fnc$distance) +#' +#' pm2 <- mergeByGroup(pm, fnc$cells) +#' +#' spe <- mergeHoodSpe(spe, pm2) +#' +#' plotHoodMat(spe, pm_cols = colnames(pm2)) +#' +setGeneric( + "plotHoodMat", + function(object, ...) standardGeneric("plotHoodMat") +) + +#' @rdname plotHoodMat +setMethod( + "plotHoodMat", + signature("matrix"), + function(object, targetCells = NA, n = 30, + hm_width = 4, hm_height = 15, clusterRows = TRUE, clusterCols = TRUE, + title = "Probability of neighborhoods") { + if (is(targetCells, "logical")) { + targetCells <- sample(rownames(object), n) + } + + dat <- object[targetCells, ] + + plotHoodMat_intl(dat, hm_width, hm_height, clusterRows, + clusterCols, title = title) + } +) + +#' @rdname plotHoodMat +setMethod( + "plotHoodMat", + signature("SpatialExperiment"), + function(object, pm_cols, targetCells = NA, n = 30, + hm_width = 4, hm_height = 15, clusterRows = TRUE, clusterCols = TRUE, + title = "Probability of neighborhoods") { + dat <- as.data.frame(colData(object)) + dat <- dat[, pm_cols] + + if (is(targetCells, "logical")) { + targetCells <- sample(rownames(dat), n) + } + + dat <- as.matrix(dat[targetCells, ]) + + plotHoodMat_intl(dat, hm_width, hm_height, clusterRows, + clusterCols, title = title) + } +) + +plotHoodMat_intl <- function(dat, hm_width, hm_height, clusterRows, + clusterCols, title = title) { + cp <- scico::scico(200, palette = "acton", direction = -1) + + ht <- ComplexHeatmap::Heatmap(dat, + col = cp, rect_gp = grid::gpar(col = "black", lwd = 1), + name = "Probability", column_title = title, + column_title_gp = grid::gpar(fontsize = 15, fontface = "italic"), + column_dend_height = unit(1, "cm"), + row_dend_width = unit(1, "cm"), column_names_rot = 45, + width = unit(hm_width, "cm"), height = unit(hm_height, "cm"), + cluster_rows = clusterRows, cluster_columns = clusterCols + ) + + return(ht) +} diff --git a/R/plot_tissue.R b/R/plot_tissue.R new file mode 100644 index 0000000..c470f4c --- /dev/null +++ b/R/plot_tissue.R @@ -0,0 +1,112 @@ +#' Plot cells based on cell position on tissue. +#' +#' @param spe SpatialExperiment object. +#' @param targetcell Optional. Can input ONE specific cell id to zoom-in on +#' the region of a specific cell. +#' @param k_near Optional. If targetcell is specified, the k_near cells around +#' the targetcell will be plotted. +#' @param targetsize Dot size of the targetcell. +#' @param targetshape Shape of the targetcell. +#' @param targetcolor Colour of the targetcell. +#' @param scaleFactor Scale factor to align with the image. +#' @param reverseY Reverse y coordinates. +#' @param ... aesthetic mappings to pass to `ggplot2::aes_string()`. +#' +#' @return A ggplot object. +#' @export +#' +#' @examples +#' +#' data("spe_test") +#' +#' plotTissue(spe, color = celltypes) +#' +plotTissue <- function(spe, targetcell = FALSE, k_near = 100, targetsize = 3, + targetshape = 1, targetcolor = "red", + scaleFactor = 1, reverseY = TRUE, ...) { + if (!is(targetcell, "logical")) { + if (length(targetcell) == 1) { + ncells <- findNearCells(spe, + k = k_near, targetCell = targetcell, reportCellID = TRUE, + reportDist = FALSE + ) |> + unlist() + + spe <- spe[, c(targetcell, ncells)] + } + } + + toplot <- SpatialExperiment::spatialCoords(spe) |> + as.data.frame() + + colnames(toplot) <- c("x", "y") + + + + cdata <- as.data.frame(SummarizedExperiment::colData(spe)) + + if ("cell_id" %in% colnames(cdata)) { + cdata <- dplyr::select(cdata, -cell_id) + } + + toplot <- toplot |> + cbind(cdata) |> + tibble::rownames_to_column("cell_id") + + aesmap <- rlang::enquos(...) + + # split aes params into those that are not aes i.e. static parametrisation + if (length(aesmap) > 0) { + is_aes <- vapply(aesmap, rlang::quo_is_symbolic, FUN.VALUE = logical(1)) + defaultmap <- lapply(aesmap[!is_aes], rlang::eval_tidy) + aesmap <- aesmap[is_aes] + } else { + defaultmap <- list() + } + + toplot[, "x"] <- scaleFactor * toplot[, "x"] + toplot[, "y"] <- scaleFactor * toplot[, "y"] + + if (reverseY) { + y_tmp <- toplot[, "y"] + mid_y <- (max(y_tmp) + min(y_tmp)) / 2 + final_y <- 2 * mid_y - y_tmp + toplot[, "y"] <- final_y + } + + p <- ggplot2::ggplot(toplot, aes(x = x, y = y, !!!aesmap)) + + p <- p + + do.call(ggplot2::geom_point, defaultmap) + + tissue_theme() + + if (!is(targetcell, "logical")) { + target_df <- toplot |> + dplyr::filter(cell_id %in% targetcell) + + p <- p + + ggplot2::geom_point(data = target_df, aes(x, y), size = targetsize, + shape = targetshape, color = targetcolor) + } + + return(p) +} + + + +tissue_theme <- function(textScale = 1.1) { + stopifnot(textScale > 0) + ggplot2::theme_minimal() + + ggplot2::theme( + panel.border = element_rect(colour = "black", fill = NA, linewidth = 1), + panel.grid = element_blank(), + legend.text = element_text(size = rel(textScale)), + legend.title = element_text(size = rel(textScale), face = "italic"), + axis.title = element_blank(), + axis.text = element_blank(), + axis.ticks = element_blank() + ) +} + + +utils::globalVariables(c("x", "y", "cell_id")) diff --git a/R/pm_clust.R b/R/pm_clust.R new file mode 100644 index 0000000..9649280 --- /dev/null +++ b/R/pm_clust.R @@ -0,0 +1,82 @@ +#' Cluster the probability matrix with K-means +#' +#' @param object A probability matrix or a SpatialExperiment. +#' @param k The number of clusters. By default is 2^ncol(object)-1. +#' @param iter.max the maximum number of iterations allowed. +#' @param nstart how many random sets should be chosen. +#' @param pm_cols The colnames of probability matrix. This is requires for +#' SpatialExperiment input. Assuming that the probability is stored in the colData. +#' @param algo Algorithm to be used. Options include Hartigan-Wong, Lloyd, and MacQueen. +#' @param ... Ignore parameter. +#' +#' @return A probability matrix or a SpatialExperiment object. For latter, +#' the clustering results are saved in the colData of the SpatialExperiment object. +#' @export +#' +#' @docType methods +#' @name clustByHood +#' @rdname clustByHood +#' @aliases clustByHood clustByHood,matrix-method clustByHood,SpatialExperiment-method +#' +#' @examples +#' +#' m <- matrix(abs(rnorm(1000 * 100)), 1000, 100) +#' +#' clust <- clustByHood(m, k = 3) +#' +setGeneric( + "clustByHood", + function(object, ...) standardGeneric("clustByHood") +) + + +#' @rdname clustByHood +setMethod( + "clustByHood", + signature("matrix"), + function(object, k = 2^ncol(object) - 1, + iter.max = 1000, nstart = 5) { + clust <- pmclust_intl(object, k, iter.max, nstart, algo = "Hartigan-Wong") + + return(clust) + } +) + +#' @rdname clustByHood +setMethod( + "clustByHood", + signature("SpatialExperiment"), + function(object, pm_cols, k = 0, + iter.max = 1000, nstart = 5, + algo = "Hartigan-Wong") { + dat <- as.data.frame(colData(object)) + + if (!all(pm_cols %in% colnames(dat))) { + stop("The pm_cols are not included in the SpatialExperiment.") + } + + dat <- dat[, pm_cols] + + if (k == 0) { + k <- 2^ncol(dat) - 1 + } + + clust <- pmclust_intl(dat, k, iter.max, nstart, algo) + + colData(object)$clusters <- as.character(clust$cluster) + metadata(object)$centroids <- clust$centers + + return(object) + } +) + + +pmclust_intl <- function(x, k, im, ns, algo) { + clust <- stats::kmeans( + x = x, centers = k, + iter.max = im, nstart = ns, + algorithm = algo + ) + + return(clust) +} diff --git a/R/read_data.R b/R/read_data.R new file mode 100644 index 0000000..7f32e85 --- /dev/null +++ b/R/read_data.R @@ -0,0 +1,106 @@ +#' Read cellular position and annotation data into a list object. +#' +#' @param spe SpatialExperiment object. +#' @param anno_col Character. The column name of the annotation to be used in +#' the following neighbourhood analysis. +#' @param cell_pos_dat data.frame object contains the cellular positions. +#' @param cell_anno_dat data.frame object contains the cell annotations. +#' @param pos_col Character. If the x and y are in the colData instead of in the +#' SpatialCoords of spe, can specify this parameter. +#' +#' @return A SpatialExperiment object. +#' +#' @export +#' +#' @examples +#' data("spe_test") +#' +#' spe <- readHoodData(spe, anno_col = "celltypes") +#' +readHoodData <- function(spe = NA, anno_col = NA, + cell_pos_dat = NA, cell_anno_dat = NA, + pos_col = NA) { + if (is(spe, "logical") & is(anno_col, "logical")) { + if (is(cell_pos_dat, "logical") | is(cell_anno_dat, "logical")) { + stop("You need to input either a SummarizedExperiment object with parameter + spe or two dataframes with parameters cell_pos_dat and cell_anno_dat") + } else if (!is.data.frame(cell_pos_dat) | !is.data.frame(cell_anno_dat)) { + stop("cell_pos_dat and cell_anno_dat should be data.frame objects.") + } + if (ncol(cell_pos_dat) != 3) { + stop("The cell_pos_dat is expected to have three columns. cell_id, x, and y.") + } + if (ncol(cell_anno_dat) != 2) { + stop("The cell_pos_dat is expected to have two columns. + cell_id and annotations.") + } + if (nrow(cell_pos_dat) != nrow(cell_anno_dat)) { + stop("The cell_pos_dat should have the same amount of cells as cell_anno_data.") + } + } else if (!is(spe, "SummarizedExperiment") | !is(anno_col, "character")) { + stop("spe should be a SummarizedExperiment object and anno_col should be + a column name in the colData of spe.") + } + + if (is(spe, "logical")) { + colnames(cell_pos_dat) <- c("cell_id", "x", "y") + colnames(cell_anno_dat) <- c("cell_id", "cell_annotation") + + dat <- list(cell_pos_dat, cell_anno_dat) + names(dat) <- c("cell_position", "cell_annotation") + + + dummy_exp_m <- matrix(1, nrow = 10, ncol = nrow(cell_pos_dat)) + colnames(dummy_exp_m) <- cell_pos_dat$cell_id + rownames(dummy_exp_m) <- paste0("gene_", seq(10)) + spe_n <- SpatialExperiment::SpatialExperiment( + assay = list("counts" = as.data.frame(dummy_exp_m)), + colData = cell_anno_dat |> tibble::column_to_rownames("cell_id"), + spatialCoords = cell_pos_dat |> tibble::column_to_rownames("cell_id") |> + as.matrix(), + metadata = list("dummy" = 1) + ) + } else { # if input is spe + if (!(anno_col %in% colnames(colData(spe)))) { + stop("anno_col is not in the colData of spe.") + } + + col_dat <- as.data.frame(colData(spe)) + colnames(col_dat)[colnames(col_dat) == anno_col] <- "cell_annotation" + + if (is(spe, "SpatialExperiment")) { + if (is(pos_col, "logical")) { + pos_dat <- spatialCoords(spe) + colnames(pos_dat) <- c("x", "y") + } else { + if (!all(pos_col %in% colnames(colData(spe)))) { + stop("pos_col is not in the colData of spe.") + } + pos_dat <- as.data.frame(colData(spe))[, pos_col] + colnames(pos_dat) <- c("x", "y") + pos_dat <- as.matrix(pos_dat) + } + } else { + if (is(pos_col, "logical") | length(pos_col) != 2) { + stop("You need to specify 2 cell position columns in colData of the spe + with parameter pos_col.") + } + + if (!all(pos_col %in% colnames(colData(spe)))) { + stop("pos_col is not in the colData of spe.") + } + pos_dat <- as.data.frame(colData(spe))[, pos_col] + colnames(pos_dat) <- c("x", "y") + pos_dat <- as.matrix(pos_dat) + } + + spe_n <- SpatialExperiment::SpatialExperiment( + assay = list("counts" = as.data.frame(as.matrix(assay(spe, "counts")))), + colData = col_dat, + spatialCoords = pos_dat, + metadata = list("dummy" = 0) + ) + } + + return(spe_n) +} diff --git a/R/soft_neighbourhood.R b/R/soft_neighbourhood.R new file mode 100644 index 0000000..128f5b6 --- /dev/null +++ b/R/soft_neighbourhood.R @@ -0,0 +1,68 @@ +#' Scan cellular neighbourhoods. +#' +#' @param m Distance matrix. Can be obtained from function findNearCells. +#' @param mode Character. Either proximityFocused or smoothFadeout. +#' By default is proximityFocused. +#' @param tau The hyperparameter tau, by default is median(m**2)/5 +#' @param t_init An initial tau. In the smoothFadeout mode, user can provide +#' an initial tau for optimization. +#' +#' @return A probability matrix. +#' @export +#' +#' @examples +#' +#' m <- matrix(abs(rnorm(1000 * 100)), 1000, 100) +#' +#' pm <- scanHoods(m) +#' +scanHoods <- function(m, mode = c("proximityFocused", "smoothFadeout"), + tau = NA, t_init = NA) { + if (length(mode) == 2) { + mode <- "proximityFocused" + } + + if (!(all(mode %in% c("proximityFocused", "smoothFadeout")))) { + stop("mode must be either proximityFocused or smoothFadeout.") + } + + if (mode == "proximityFocused") { + if (is(tau, "logical")) { + tau <- median(m**2) / 5 + } + } else if (mode == "smoothFadeout") { + if (is(t_init, "logical")) { + t_init <- median(m**2) + } + + optim_result <- stats::optim(par = t_init, fn = f_nll, + m = m, method = "BFGS") + + tau <- optim_result$par + + msg <- paste0("Optimized tau is: ", tau) + + message(msg) + } + + msg <- paste0("Tau is set to: ", tau) + message(msg) + probm <- soft_max_intl(m, tau) + + return(probm) +} + + +soft_max_intl <- function(m, t) { + m <- m**2 + m <- -m / t + exp_m <- exp(m) + probm <- sweep(exp_m, 1, rowSums(exp_m), "/") + return(probm) +} + +f_nll <- function(m, t) { + probm <- soft_max_intl(m, t) + nll <- -sum(log(probm + 1e-8)) + return(nll) +} diff --git a/README.md b/README.md new file mode 100644 index 0000000..b773175 --- /dev/null +++ b/README.md @@ -0,0 +1,16 @@ +# hoodscanR: Cellular neighborhoods scanning from spatial transcriptomics data in R + +[![R-CMD-check](https://github.com/DavisLaboratory/hoodscanR/workflows/R-CMD-check-bioc/badge.svg)](https://github.com/DavisLaboratory/hoodscanR/actions) + + + + +Install development version from GitHub + +``` +library(devtools) +devtools::install_github("DavisLaboratory/hoodscanR") +``` + + + diff --git a/data/spe_test.rda b/data/spe_test.rda new file mode 100644 index 0000000..285e25c Binary files /dev/null and b/data/spe_test.rda differ diff --git a/man/calcMetrics.Rd b/man/calcMetrics.Rd new file mode 100644 index 0000000..2115685 --- /dev/null +++ b/man/calcMetrics.Rd @@ -0,0 +1,40 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/calc_metrics.R +\name{calcMetrics} +\alias{calcMetrics} +\title{Calculate metrics for probability matrix} +\usage{ +calcMetrics(spe, pm = NA, pm_cols = NA) +} +\arguments{ +\item{spe}{A SpatialExperiment object.} + +\item{pm}{Optional. The probability matrix.} + +\item{pm_cols}{The colnames of probability matrix. This is requires for +SpatialExperiment input. Assuming that the probability is stored in the colData.} +} +\value{ +A SpatialExperiment object. Calculated entropy and perplexity are +saved as columns in the colData of the SpatialExperiment object. +} +\description{ +Calculate metrics for probability matrix +} +\examples{ + +data("spe_test") + +spe <- readHoodData(spe, anno_col = "celltypes") + +fnc <- findNearCells(spe, k = 100) + +pm <- scanHoods(fnc$distance) + +pm2 <- mergeByGroup(pm, fnc$cells) + +spe <- mergeHoodSpe(spe, pm2) + +spe <- calcMetrics(spe, pm_cols = colnames(pm2)) + +} diff --git a/man/clustByHood.Rd b/man/clustByHood.Rd new file mode 100644 index 0000000..dd0fa3f --- /dev/null +++ b/man/clustByHood.Rd @@ -0,0 +1,52 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/pm_clust.R +\docType{methods} +\name{clustByHood} +\alias{clustByHood} +\alias{clustByHood,matrix-method} +\alias{clustByHood,SpatialExperiment-method} +\title{Cluster the probability matrix with K-means} +\usage{ +clustByHood(object, ...) + +\S4method{clustByHood}{matrix}(object, k = 2^ncol(object) - 1, iter.max = 1000, nstart = 5) + +\S4method{clustByHood}{SpatialExperiment}( + object, + pm_cols, + k = 0, + iter.max = 1000, + nstart = 5, + algo = "Hartigan-Wong" +) +} +\arguments{ +\item{object}{A probability matrix or a SpatialExperiment.} + +\item{...}{Ignore parameter.} + +\item{k}{The number of clusters. By default is 2^ncol(object)-1.} + +\item{iter.max}{the maximum number of iterations allowed.} + +\item{nstart}{how many random sets should be chosen.} + +\item{pm_cols}{The colnames of probability matrix. This is requires for +SpatialExperiment input. Assuming that the probability is stored in the colData.} + +\item{algo}{Algorithm to be used. Options include Hartigan-Wong, Lloyd, and MacQueen.} +} +\value{ +A probability matrix or a SpatialExperiment object. For latter, +the clustering results are saved in the colData of the SpatialExperiment object. +} +\description{ +Cluster the probability matrix with K-means +} +\examples{ + +m <- matrix(abs(rnorm(1000 * 100)), 1000, 100) + +clust <- clustByHood(m, k = 3) + +} diff --git a/man/figures/hoodscanR_sticker.png b/man/figures/hoodscanR_sticker.png new file mode 100644 index 0000000..e653f3f Binary files /dev/null and b/man/figures/hoodscanR_sticker.png differ diff --git a/man/findNearCells.Rd b/man/findNearCells.Rd new file mode 100644 index 0000000..bff70d3 --- /dev/null +++ b/man/findNearCells.Rd @@ -0,0 +1,49 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/find_near_cells.R +\name{findNearCells} +\alias{findNearCells} +\title{Find the k-th nearest cells for each cell} +\usage{ +findNearCells( + dat, + k = 100, + targetCell = FALSE, + reportCellID = FALSE, + reportDist = TRUE, + anno_col = 0 +) +} +\arguments{ +\item{dat}{A SpatialExperiment object, can be generated using function \code{readHoodData}.} + +\item{k}{The maximum number of nearest cells to compute.} + +\item{targetCell}{Specify the cells to be the target cell for finding nearest cells.} + +\item{reportCellID}{Logical. Set to TRUE to report cell id instead of cell types.} + +\item{reportDist}{Logical. Set to TRUE to report the distance matrix.} + +\item{anno_col}{Character vector. The name of annotation column to use.} +} +\value{ +A list includes a data.frame and a matrix, describing the cell types +and distances of the k-th nearest cells of each cell. +} +\description{ +Find the k-th nearest cells for each cell +} +\details{ +The \code{findNearCells} function uses the \code{nn2} function from the \code{RANN} package, +which uses the Approximate Near Neighbor (ANN) C++ library. For more infromation on the +ANN library please see http://www.cs.umd.edu/~mount/ANN/. +} +\examples{ + +data("spe_test") + +spe <- readHoodData(spe, anno_col = "celltypes") + +fnc <- findNearCells(spe, k = 100) + +} diff --git a/man/hoodscanR-package.Rd b/man/hoodscanR-package.Rd new file mode 100644 index 0000000..f112b45 --- /dev/null +++ b/man/hoodscanR-package.Rd @@ -0,0 +1,18 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/hoodscanR-package.R +\docType{package} +\name{hoodscanR-package} +\alias{hoodscanR-package} +\alias{hoodscanR} +\title{Method to identify cellular spatial neighbourhood from single cell spatial +transcriptomics data.} +\description{ +\code{hoodscanR} implements a novel method to scan for cell neighbourhood from +spatial transcriptomics data at single cell level, such as CosMx and MERFISH etc. +\code{hoodscanR} takes the cellular position and cell type annotations +as inputs, allowing cellular spatial neighbourhood analysis. +} +\author{ +Ning Liu \email{liu.n@wehi.edu.au} +} +\keyword{internal} diff --git a/man/mergeByGroup.Rd b/man/mergeByGroup.Rd new file mode 100644 index 0000000..5d86576 --- /dev/null +++ b/man/mergeByGroup.Rd @@ -0,0 +1,33 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/merge_by_group.R +\name{mergeByGroup} +\alias{mergeByGroup} +\title{Merge probability matrix based on annotations} +\usage{ +mergeByGroup(pm, group_df) +} +\arguments{ +\item{pm}{A numeric matrix. Probability matrix generated by the soft_max function.} + +\item{group_df}{A character matrix. Annotation of the neighboring cells to be used.} +} +\value{ +A probability matrix, describing the probability of each cell being +in each cellular neighborhood. +} +\description{ +Merge probability matrix based on annotations +} +\examples{ + +data("spe_test") + +spe <- readHoodData(spe, anno_col = "celltypes") + +fnc <- findNearCells(spe, k = 100) + +pm <- scanHoods(fnc$distance) + +pm2 <- mergeByGroup(pm, fnc$cells) + +} diff --git a/man/mergeHoodSpe.Rd b/man/mergeHoodSpe.Rd new file mode 100644 index 0000000..ced232f --- /dev/null +++ b/man/mergeHoodSpe.Rd @@ -0,0 +1,35 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/merge_pm_spe.R +\name{mergeHoodSpe} +\alias{mergeHoodSpe} +\title{Merge probability matrix into SpatialExperiment object.} +\usage{ +mergeHoodSpe(spe, pm) +} +\arguments{ +\item{spe}{A SpatialExperiment object.} + +\item{pm}{Probability matrix. Can be obtained by the function mergeByGroup.} +} +\value{ +A SpatialExperiment object. Cell-level neighborhood information are +saved in the colData of the SpatialExperiment object. +} +\description{ +Merge probability matrix into SpatialExperiment object. +} +\examples{ + +data("spe_test") + +spe <- readHoodData(spe, anno_col = "celltypes") + +fnc <- findNearCells(spe, k = 100) + +pm <- scanHoods(fnc$distance) + +pm2 <- mergeByGroup(pm, fnc$cells) + +spe <- mergeHoodSpe(spe, pm2) + +} diff --git a/man/plotColocal.Rd b/man/plotColocal.Rd new file mode 100644 index 0000000..f8f1999 --- /dev/null +++ b/man/plotColocal.Rd @@ -0,0 +1,77 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/plot_heatmap.R +\docType{methods} +\name{plotColocal} +\alias{plotColocal} +\alias{plotColocal,} +\alias{matrix-method} +\alias{SpatialExperiment-method} +\alias{plotColocal,matrix-method} +\alias{plotColocal,SpatialExperiment-method} +\title{Plot heatmap for neighbourhood analysis} +\usage{ +plotColocal(object, ...) + +\S4method{plotColocal}{matrix}(object, hm_width = 5, hm_height = 5) + +\S4method{plotColocal}{SpatialExperiment}( + object, + pm_cols, + self.cor = TRUE, + byGroup = NA, + hm_width = 5, + hm_height = 5, + clusterRows = TRUE, + clusterCols = TRUE, + return_matrix = FALSE +) +} +\arguments{ +\item{object}{A probability matrix or SpatialExperiment.} + +\item{...}{Ignore parameter.} + +\item{hm_width}{Integer. The width of heatmap.} + +\item{hm_height}{Integer. The height of heatmap.} + +\item{pm_cols}{The colnames of probability matrix. This is requires for +SpatialExperiment input. Assuming that the probability is stored in the colData.} + +\item{self.cor}{Logical. By default is TRUE, inidicating running a correlation +between neighbourhoods to perform a simple co-localization analysis. +When this set to FALSE, it will plot the average probability of each +neighbourhood by group using the byGroup parameter.} + +\item{byGroup}{Character. This is required when self.cor is set to FALSE.} + +\item{clusterRows}{Logical. Cluster rows.} + +\item{clusterCols}{Logical. Cluster columns.} + +\item{return_matrix}{Logical. Export a numeric matrix .} +} +\value{ +A ComplexHeatmap plot. When return_matrix is set to TRUE, +return a matrix Object. +} +\description{ +Plot heatmap for neighbourhood analysis +} +\examples{ + +data("spe_test") + +spe <- readHoodData(spe, anno_col = "celltypes") + +fnc <- findNearCells(spe, k = 100) + +pm <- scanHoods(fnc$distance) + +pm2 <- mergeByGroup(pm, fnc$cells) + +spe <- mergeHoodSpe(spe, pm2) + +plotColocal(spe, pm_cols = colnames(pm2)) + +} diff --git a/man/plotHoodMat.Rd b/man/plotHoodMat.Rd new file mode 100644 index 0000000..253c7c0 --- /dev/null +++ b/man/plotHoodMat.Rd @@ -0,0 +1,80 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/plot_pm.R +\docType{methods} +\name{plotHoodMat} +\alias{plotHoodMat} +\alias{plotHoodMat,matrix-method} +\alias{plotHoodMat,SpatialExperiment-method} +\title{Plot probability matrix as a heatmap} +\usage{ +plotHoodMat(object, ...) + +\S4method{plotHoodMat}{matrix}( + object, + targetCells = NA, + n = 30, + hm_width = 4, + hm_height = 15, + clusterRows = TRUE, + clusterCols = TRUE, + title = "Probability of neighborhoods" +) + +\S4method{plotHoodMat}{SpatialExperiment}( + object, + pm_cols, + targetCells = NA, + n = 30, + hm_width = 4, + hm_height = 15, + clusterRows = TRUE, + clusterCols = TRUE, + title = "Probability of neighborhoods" +) +} +\arguments{ +\item{object}{A probability matrix or SpatialExperiment.} + +\item{...}{Ignore parameter.} + +\item{targetCells}{Character. Optional. Can speicify one or more cells to be plotted.} + +\item{n}{Integer. The number of randomly selected cells to be plotted. +This parameter will be used when targetCells is not specify.} + +\item{hm_width}{Integer. The width of heatmap.} + +\item{hm_height}{Integer. The height of heatmap.} + +\item{clusterRows}{Logical. Cluster rows or not.} + +\item{clusterCols}{Logical. Cluster columns or not.} + +\item{title}{Title of the heatmap.} + +\item{pm_cols}{The colnames of probability matrix. This is requires for +SpatialExperiment input. Assuming that the probability is stored in the colData.} +} +\value{ +A ComplexHeatmap plot. +} +\description{ +Plot probability matrix as a heatmap +} +\examples{ + +data("spe_test") + +spe <- readHoodData(spe, anno_col = "celltypes") + +fnc <- findNearCells(spe, k = 100) + +pm <- scanHoods(fnc$distance) + +pm2 <- mergeByGroup(pm, fnc$cells) + +spe <- mergeHoodSpe(spe, pm2) + +plotHoodMat(spe, pm_cols = colnames(pm2)) + +} diff --git a/man/plotProbDist.Rd b/man/plotProbDist.Rd new file mode 100644 index 0000000..05b7402 --- /dev/null +++ b/man/plotProbDist.Rd @@ -0,0 +1,67 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/plot_pd.R +\docType{methods} +\name{plotProbDist} +\alias{plotProbDist} +\alias{plotProbDist,matrix-method} +\alias{plotProbDist,SpatialExperiment-method} +\title{Plot probability distribution} +\usage{ +plotProbDist(object, ...) + +\S4method{plotProbDist}{matrix}(object, targetCells = NA, ...) + +\S4method{plotProbDist}{SpatialExperiment}( + object, + pm_cols, + targetCells = NA, + by_cluster = FALSE, + show_clusters = as.character(seq(6)), + plot_all = FALSE, + sample_size = 2, + ... +) +} +\arguments{ +\item{object}{A probability matrix or SpatialExperiment.} + +\item{...}{aesthetic mappings to pass to \code{ggplot2::aes_string()}.} + +\item{targetCells}{Character. Optional. Can speicify one or more cells to be plotted.} + +\item{pm_cols}{The colnames of probability matrix. This is requires for +SpatialExperiment input. Assuming that the probability is stored in the colData.} + +\item{by_cluster}{Logical. By default is TRUE, to plot distribution by each cluster.} + +\item{show_clusters}{Character. The cluster to be ploted, by default is 1 to 6.} + +\item{plot_all}{Logical. By default is FALSE, set this to true to plot box +plot instead of bar plot to show all cells in each cluster.} + +\item{sample_size}{Integer. By default is 2, sampling two cell from each +cluster to be plotted.} +} +\value{ +A ggplot object. +} +\description{ +Plot probability distribution +} +\examples{ + +data("spe_test") + +spe <- readHoodData(spe, anno_col = "celltypes") + +fnc <- findNearCells(spe, k = 100) + +pm <- scanHoods(fnc$distance) + +pm2 <- mergeByGroup(pm, fnc$cells) + +spe <- mergeHoodSpe(spe, pm2) + +plotProbDist(spe, pm_cols = colnames(pm2)) + +} diff --git a/man/plotTissue.Rd b/man/plotTissue.Rd new file mode 100644 index 0000000..7ebbc60 --- /dev/null +++ b/man/plotTissue.Rd @@ -0,0 +1,52 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/plot_tissue.R +\name{plotTissue} +\alias{plotTissue} +\title{Plot cells based on cell position on tissue.} +\usage{ +plotTissue( + spe, + targetcell = FALSE, + k_near = 100, + targetsize = 3, + targetshape = 1, + targetcolor = "red", + scaleFactor = 1, + reverseY = TRUE, + ... +) +} +\arguments{ +\item{spe}{SpatialExperiment object.} + +\item{targetcell}{Optional. Can input ONE specific cell id to zoom-in on +the region of a specific cell.} + +\item{k_near}{Optional. If targetcell is specified, the k_near cells around +the targetcell will be plotted.} + +\item{targetsize}{Dot size of the targetcell.} + +\item{targetshape}{Shape of the targetcell.} + +\item{targetcolor}{Colour of the targetcell.} + +\item{scaleFactor}{Scale factor to align with the image.} + +\item{reverseY}{Reverse y coordinates.} + +\item{...}{aesthetic mappings to pass to \code{ggplot2::aes_string()}.} +} +\value{ +A ggplot object. +} +\description{ +Plot cells based on cell position on tissue. +} +\examples{ + +data("spe_test") + +plotTissue(spe, color = celltypes) + +} diff --git a/man/readHoodData.Rd b/man/readHoodData.Rd new file mode 100644 index 0000000..7e45248 --- /dev/null +++ b/man/readHoodData.Rd @@ -0,0 +1,39 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/read_data.R +\name{readHoodData} +\alias{readHoodData} +\title{Read cellular position and annotation data into a list object.} +\usage{ +readHoodData( + spe = NA, + anno_col = NA, + cell_pos_dat = NA, + cell_anno_dat = NA, + pos_col = NA +) +} +\arguments{ +\item{spe}{SpatialExperiment object.} + +\item{anno_col}{Character. The column name of the annotation to be used in +the following neighbourhood analysis.} + +\item{cell_pos_dat}{data.frame object contains the cellular positions.} + +\item{cell_anno_dat}{data.frame object contains the cell annotations.} + +\item{pos_col}{Character. If the x and y are in the colData instead of in the +SpatialCoords of spe, can specify this parameter.} +} +\value{ +A SpatialExperiment object. +} +\description{ +Read cellular position and annotation data into a list object. +} +\examples{ +data("spe_test") + +spe <- readHoodData(spe, anno_col = "celltypes") + +} diff --git a/man/scanHoods.Rd b/man/scanHoods.Rd new file mode 100644 index 0000000..38da2df --- /dev/null +++ b/man/scanHoods.Rd @@ -0,0 +1,37 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/soft_neighbourhood.R +\name{scanHoods} +\alias{scanHoods} +\title{Scan cellular neighbourhoods.} +\usage{ +scanHoods( + m, + mode = c("proximityFocused", "smoothFadeout"), + tau = NA, + t_init = NA +) +} +\arguments{ +\item{m}{Distance matrix. Can be obtained from function findNearCells.} + +\item{mode}{Character. Either proximityFocused or smoothFadeout. +By default is proximityFocused.} + +\item{tau}{The hyperparameter tau, by default is median(m**2)/5} + +\item{t_init}{An initial tau. In the smoothFadeout mode, user can provide +an initial tau for optimization.} +} +\value{ +A probability matrix. +} +\description{ +Scan cellular neighbourhoods. +} +\examples{ + +m <- matrix(abs(rnorm(1000 * 100)), 1000, 100) + +pm <- scanHoods(m) + +} diff --git a/man/spe_test.Rd b/man/spe_test.Rd new file mode 100644 index 0000000..99c07fa --- /dev/null +++ b/man/spe_test.Rd @@ -0,0 +1,25 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/hoodscanR-package.R +\docType{data} +\name{spe_test} +\alias{spe_test} +\alias{spe} +\title{Example test data} +\format{ +A SpatialExperiment object +} +\usage{ +data("spe_test") +} +\value{ +A SpatialExperiment object +} +\description{ +hoodscanR-package has 1 datasets: \itemize{ +\item spe_test An example data +} +} +\examples{ +data(spe_test) +} +\keyword{internal} diff --git a/src/RcppExports.cpp b/src/RcppExports.cpp new file mode 100644 index 0000000..a3d4d95 --- /dev/null +++ b/src/RcppExports.cpp @@ -0,0 +1,33 @@ +// Generated by using Rcpp::compileAttributes() -> do not edit by hand +// Generator token: 10BE3573-1514-4C36-9D1C-5A225CD40393 + +#include + +using namespace Rcpp; + +#ifdef RCPP_USE_GLOBAL_ROSTREAM +Rcpp::Rostream& Rcpp::Rcout = Rcpp::Rcpp_cout_get(); +Rcpp::Rostream& Rcpp::Rcerr = Rcpp::Rcpp_cerr_get(); +#endif + +// calculate_metrics +DataFrame calculate_metrics(NumericMatrix probs); +RcppExport SEXP _hoodscanR_calculate_metrics(SEXP probsSEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< NumericMatrix >::type probs(probsSEXP); + rcpp_result_gen = Rcpp::wrap(calculate_metrics(probs)); + return rcpp_result_gen; +END_RCPP +} + +static const R_CallMethodDef CallEntries[] = { + {"_hoodscanR_calculate_metrics", (DL_FUNC) &_hoodscanR_calculate_metrics, 1}, + {NULL, NULL, 0} +}; + +RcppExport void R_init_hoodscanR(DllInfo *dll) { + R_registerRoutines(dll, NULL, CallEntries, NULL, NULL); + R_useDynamicSymbols(dll, FALSE); +} diff --git a/src/RcppExports.o b/src/RcppExports.o new file mode 100644 index 0000000..7f58ad6 Binary files /dev/null and b/src/RcppExports.o differ diff --git a/src/bhattacharyya_distance.o b/src/bhattacharyya_distance.o new file mode 100644 index 0000000..b275b19 Binary files /dev/null and b/src/bhattacharyya_distance.o differ diff --git a/src/cal_metrics.cpp b/src/cal_metrics.cpp new file mode 100644 index 0000000..e205d53 --- /dev/null +++ b/src/cal_metrics.cpp @@ -0,0 +1,57 @@ +#include +#include + +using namespace Rcpp; + +// [[Rcpp::export]] +DataFrame calculate_metrics(NumericMatrix probs) { + + int n_rows = probs.nrow(); + int n_cols = probs.ncol(); + + NumericVector entropy(n_rows); + NumericVector perplexity(n_rows); + //NumericVector mean(n_rows); + //NumericVector variance(n_rows); + //NumericVector sd(n_rows); + + for (int i = 0; i < n_rows; i++) { + + // calculate entropy + entropy[i] = 0; + for (int j = 0; j < n_cols; j++) { + if (probs(i, j) > 0) { + entropy[i] -= probs(i, j) * log2(probs(i, j)); + } + } + + // calculate perplexity + perplexity[i] = pow(2, entropy[i]); + + // calculate mean + //mean[i] = 0; + //for (int j = 0; j < n_cols; j++) { + // mean[i] += (j+1) * probs(i, j); + //} + + // calculate variance + //variance[i] = 0; + //for (int j = 0; j < n_cols; j++) { + // variance[i] += pow((j+1) - mean[i], 2) * probs(i, j); + //} + + // calculate standard deviation + //sd[i] = sqrt(variance[i]); + + } + // create output data frame + DataFrame output = DataFrame::create( + _["entropy"] = entropy, + _["perplexity"] = perplexity + //_["mean"] = mean, + //_["variance"] = variance, + //_["sd"] = sd + ); + + return output; +} diff --git a/src/cal_metrics.o b/src/cal_metrics.o new file mode 100644 index 0000000..9725861 Binary files /dev/null and b/src/cal_metrics.o differ diff --git a/src/extra_norm.o b/src/extra_norm.o new file mode 100644 index 0000000..b3238bd Binary files /dev/null and b/src/extra_norm.o differ diff --git a/src/hoodscanR.so b/src/hoodscanR.so new file mode 100755 index 0000000..7927878 Binary files /dev/null and b/src/hoodscanR.so differ diff --git a/src/kl_divergence.o b/src/kl_divergence.o new file mode 100644 index 0000000..cb6105f Binary files /dev/null and b/src/kl_divergence.o differ diff --git a/src/soft_max_intl.o b/src/soft_max_intl.o new file mode 100644 index 0000000..89da5d4 Binary files /dev/null and b/src/soft_max_intl.o differ diff --git a/src/softmax_loop.o b/src/softmax_loop.o new file mode 100644 index 0000000..3fe46c1 Binary files /dev/null and b/src/softmax_loop.o differ diff --git a/tests/testthat.R b/tests/testthat.R new file mode 100644 index 0000000..ce79249 --- /dev/null +++ b/tests/testthat.R @@ -0,0 +1,12 @@ +# This file is part of the standard setup for testthat. +# It is recommended that you do not modify it. +# +# Where should you do additional test configuration? +# Learn more about the roles of various files in: +# * https://r-pkgs.org/tests.html +# * https://testthat.r-lib.org/reference/test_package.html#special-files + +library(testthat) +library(hoodscanR) + +test_check("hoodscanR") diff --git a/tests/testthat/test-calcMetrics.R b/tests/testthat/test-calcMetrics.R new file mode 100644 index 0000000..2f14eb6 --- /dev/null +++ b/tests/testthat/test-calcMetrics.R @@ -0,0 +1,16 @@ +test_that("Testing calcMetrics", { + data("spe_test") + spe <- readHoodData(spe, anno_col = "celltypes") + fnc <- findNearCells(spe, k = 100) + pm <- scanHoods(fnc$distance) + pm2 <- mergeByGroup(pm, fnc$cells) + spe <- mergeHoodSpe(spe, pm2) + + expect_false("entropy" %in% colnames(colData(spe))) + expect_false("perplexity" %in% colnames(colData(spe))) + + spe <- calcMetrics(spe, pm_cols = colnames(pm2)) + + expect_true("entropy" %in% colnames(colData(spe))) + expect_true("perplexity" %in% colnames(colData(spe))) +}) diff --git a/tests/testthat/test-clustByHood.R b/tests/testthat/test-clustByHood.R new file mode 100644 index 0000000..6da3669 --- /dev/null +++ b/tests/testthat/test-clustByHood.R @@ -0,0 +1,20 @@ +test_that("Testing clustByHood", { + m <- matrix(abs(rnorm(1000 * 100)), 1000, 100) + + clust <- clustByHood(m, k = 3) + + expect_equal(class(clust), "kmeans") + + data("spe_test") + spe <- readHoodData(spe, anno_col = "celltypes") + fnc <- findNearCells(spe, k = 100) + pm <- scanHoods(fnc$distance) + pm2 <- mergeByGroup(pm, fnc$cells) + spe <- mergeHoodSpe(spe, pm2) + + spe <- clustByHood(spe, pm_cols = colnames(pm2), k = 5) + + expect_true("clusters" %in% colnames(colData(spe))) + + expect_equal(typeof(colData(spe)$clusters), "character") +}) diff --git a/tests/testthat/test-findNearCells.R b/tests/testthat/test-findNearCells.R new file mode 100644 index 0000000..01b26fa --- /dev/null +++ b/tests/testthat/test-findNearCells.R @@ -0,0 +1,17 @@ +test_that("Testing findNearCells", { + data("spe_test") + + spe <- readHoodData(spe, anno_col = "celltypes") + + fnc <- findNearCells(spe, k = 100) + + expect_true(is.list(fnc)) + + expect_equal(length(fnc), 2) + + expect_equal(dim(fnc[[1]]), dim(fnc[[2]])) + + fnc2 <- findNearCells(spe, k = 100, reportDist = FALSE) + + expect_true(is.data.frame(fnc2)) +}) diff --git a/tests/testthat/test-mergeByGroup.R b/tests/testthat/test-mergeByGroup.R new file mode 100644 index 0000000..f356e4b --- /dev/null +++ b/tests/testthat/test-mergeByGroup.R @@ -0,0 +1,11 @@ +test_that("Testing mergeByGroup", { + data("spe_test") + spe <- readHoodData(spe, anno_col = "celltypes") + fnc <- findNearCells(spe, k = 100) + pm <- scanHoods(fnc$distance) + pm2 <- mergeByGroup(pm, fnc$cells) + + expect_equal(ncol(pm2), 6) + expect_equal(nrow(pm2), nrow(pm)) + expect_equal(as.vector(rowSums(pm2)), rep(1, nrow(pm))) +}) diff --git a/tests/testthat/test-mergeHoodSpe.R b/tests/testthat/test-mergeHoodSpe.R new file mode 100644 index 0000000..e541ba8 --- /dev/null +++ b/tests/testthat/test-mergeHoodSpe.R @@ -0,0 +1,13 @@ +test_that("Testing mergeHoodSpe", { + data("spe_test") + spe <- readHoodData(spe, anno_col = "celltypes") + fnc <- findNearCells(spe, k = 100) + pm <- scanHoods(fnc$distance) + pm2 <- mergeByGroup(pm, fnc$cells) + + expect_false(all(colnames(pm2) %in% colnames(colData(spe)))) + + spe <- mergeHoodSpe(spe, pm2) + + expect_true(all(colnames(pm2) %in% colnames(colData(spe)))) +}) diff --git a/tests/testthat/test-plotColocal.R b/tests/testthat/test-plotColocal.R new file mode 100644 index 0000000..cdc451e --- /dev/null +++ b/tests/testthat/test-plotColocal.R @@ -0,0 +1,11 @@ +test_that("Testing plotColocal", { + data("spe_test") + spe <- readHoodData(spe, anno_col = "celltypes") + fnc <- findNearCells(spe, k = 100) + pm <- scanHoods(fnc$distance) + pm2 <- mergeByGroup(pm, fnc$cells) + spe <- mergeHoodSpe(spe, pm2) + + expect_silent(plotColocal(spe, pm_cols = colnames(pm2))) + +}) diff --git a/tests/testthat/test-plotHoodMat.R b/tests/testthat/test-plotHoodMat.R new file mode 100644 index 0000000..2bddaf5 --- /dev/null +++ b/tests/testthat/test-plotHoodMat.R @@ -0,0 +1,10 @@ +test_that("testing plotHoodMat", { + data("spe_test") + spe <- readHoodData(spe, anno_col = "celltypes") + fnc <- findNearCells(spe, k = 100) + pm <- scanHoods(fnc$distance) + pm2 <- mergeByGroup(pm, fnc$cells) + spe <- mergeHoodSpe(spe, pm2) + + expect_silent(plotHoodMat(spe, pm_cols = colnames(pm2))) +}) diff --git a/tests/testthat/test-plotProbDist.R b/tests/testthat/test-plotProbDist.R new file mode 100644 index 0000000..7602398 --- /dev/null +++ b/tests/testthat/test-plotProbDist.R @@ -0,0 +1,14 @@ +test_that("Testing plotProbDist", { + data("spe_test") + spe <- readHoodData(spe, anno_col = "celltypes") + fnc <- findNearCells(spe, k = 100) + pm <- scanHoods(fnc$distance) + pm2 <- mergeByGroup(pm, fnc$cells) + spe <- mergeHoodSpe(spe, pm2) + + expect_silent(plotProbDist(spe, pm_cols = colnames(pm2))) + + spe <- clustByHood(spe, pm_cols = colnames(pm2), k = 5) + + expect_silent(plotProbDist(spe, pm_cols = colnames(pm2), by_cluster = TRUE, plot_all = TRUE)) +}) diff --git a/tests/testthat/test-plotTissue.R b/tests/testthat/test-plotTissue.R new file mode 100644 index 0000000..6111bc5 --- /dev/null +++ b/tests/testthat/test-plotTissue.R @@ -0,0 +1,5 @@ +test_that("testing plotTissue", { + data("spe_test") + + expect_silent(plotTissue(spe, color = celltypes)) +}) diff --git a/tests/testthat/test-readHoodData.R b/tests/testthat/test-readHoodData.R new file mode 100644 index 0000000..2307481 --- /dev/null +++ b/tests/testthat/test-readHoodData.R @@ -0,0 +1,10 @@ +test_that("Testing readHoodData", { + data("spe_test") + + spe1 <- readHoodData(spe, anno_col = "celltypes") + + expect_true("cell_annotation" %in% colnames(colData(spe1))) + + expect_error(readHoodData(spe, anno_col = "xxx")) + +}) diff --git a/tests/testthat/test-scanHoods.R b/tests/testthat/test-scanHoods.R new file mode 100644 index 0000000..1e688b0 --- /dev/null +++ b/tests/testthat/test-scanHoods.R @@ -0,0 +1,11 @@ +test_that("Testing scanHoods", { + m <- matrix(abs(rnorm(1000 * 100)), 1000, 100) + + pm <- scanHoods(m) + + expect_equal(nrow(pm),1000) + + expect_equal(ncol(pm), 100) + + expect_equal(rowSums(pm), rep(1, 1000)) +}) diff --git a/vignettes/Quick_start.Rmd b/vignettes/Quick_start.Rmd new file mode 100644 index 0000000..69f97ad --- /dev/null +++ b/vignettes/Quick_start.Rmd @@ -0,0 +1,186 @@ +--- +title: "A quick start guide to the standR package" +author: "Ning Liu, Melissa Davis" +date: "`r Sys.Date()`" +output: + html_document: + toc: true + number_sections: true + theme: cosmo + highlight: tango + code_folding: show +vignette: > + %\VignetteIndexEntry{hoodscanR_introduction} + %\VignetteEngine{knitr::rmarkdown} + %\VignetteEncoding{UTF-8} +--- + + + +```{r message=FALSE, warning=FALSE, include=FALSE} +knitr::opts_chunk$set(echo = TRUE, message = FALSE, warning = FALSE, dpi = 72, fig.retina = 1) +suppressWarnings(library(ggplot2)) + + +``` + +# Installation + +```{r eval=FALSE} +if (!require("BiocManager", quietly = TRUE)) { + install.packages("BiocManager") +} + +BiocManager::install("hoodscanR") +``` + +The development version of `standR` can be installed from GitHub: + +```{r eval=FALSE} +devtools::install_github("DavisLaboratory/hoodscanR") +``` + +# Quick start + +```{r} +library(hoodscanR) +library(SpatialExperiment) +library(scico) +``` + +# Data exploration + +```{r} +data("spe_test") + +spe <- readHoodData(spe, anno_col = "celltypes") +``` + +```{r} +spe +``` + +```{r} +colData(spe) +``` + +We can have a look at the tissue and cell positions by using the function `plotTissue`. + +The test data is relatively sparse with low-level cell type annotations. + +```{r} +col.pal <- c("red3", "royalblue", "gold", "cyan2", "purple3", "darkgreen") + +plotTissue(spe, color = cell_annotation, size = 1.5, alpha = 0.8) + + scale_color_manual(values = col.pal) +``` + +# Neighborhoods scanning + +In order to perform neighborhood scanning, we need to firstly identify k (in this example, k = 100) nearest cells for each cells. + +```{r} +fnc <- findNearCells(spe, k = 100) +``` + +The output of `findNearCells` function includes two matrix, an annotation matrix and a distance matrix. + +```{r} +lapply(fnc, function(x) x[1:10, 1:5]) +``` + +We can then perform neighborhood analysis using the function `scanHoods`. + +```{r} +pm <- scanHoods(fnc$distance) +``` + +The resulting matrix is the probability of each cell associating with their 100 nearest cells. + +```{r} +pm[1:10, 1:5] +``` + +We can then merge the probabilities by the cell types of the 100 nearest cells. + +```{r} +hoods <- mergeByGroup(pm, fnc$cells) +``` + +```{r} +hoods[1:10, ] +``` + +# Neighborhoods analysis + +We plot randomly plot 10 cells to see the output of neighborhood scanning using `plotHoodMat`. + +```{r} +plotHoodMat(hoods, n = 10, hm_height = 5) +``` + +Or to check the cells-of-interest with the parameter `targetCells` within the function + +```{r} +plotHoodMat(hoods, targetCells = c("Lung9_Rep1_5_1975", "Lung9_Rep1_5_2712"), hm_height = 3) +``` + +We can then merge the neighborhood results with the SpatialExperiment object using `mergeHoodSpe` so that we can conduct more neighborhood-related analysis. + +```{r} +spe <- mergeHoodSpe(spe, hoods) +``` + +To summarise our neighborhood results, we can use `calcMetrics` to calculate entropy and perplexity of the probability matrix so that we can have a summarisatino of the neighborhood distribution across the tissue slide, i.e. where neighborhood is more distinct and where is more mixed. + +```{r} +spe <- calcMetrics(spe, pm_cols = colnames(hoods)) +``` + +We then again use `plotTissue` to plot out the entropy or perplexity. + +```{r} +plotTissue(spe, size = 1.5, color = entropy) + + scale_color_scico(palette = "tokyo") +``` + +```{r} +plotTissue(spe, size = 1.5, color = perplexity) + + scale_color_scico(palette = "tokyo") +``` + +We can perform neighborhood colocalization analysis using `plotColocal`. Here we can see in the test data, endothelial cells and stromal cells are more likely to colocalize, epithelial cells and dividing cells are more likely to colocalize. + +```{r} +plotColocal(spe, pm_cols = colnames(hoods)) +``` + +We can cluster the cells by their neighborhood probability distribution using `clustByHood`, it is based on the k-means algorithm and here we set k to 10. + +```{r} +spe <- clustByHood(spe, pm_cols = colnames(hoods), k = 10) +``` + +We can see what are the neighborhood distributions look like in each cluster using `plotProbDist`. + +```{r} +plotProbDist(spe, pm_cols = colnames(hoods), by_cluster = TRUE, plot_all = TRUE, show_clusters = as.character(seq(10))) +``` + +We can plot the clusters on the tissue slide, agian using `plotTissue`. + +```{r} +plotTissue(spe, color = clusters) +``` + +# Session info + +```{r} +sessionInfo() +``` + +