### Install bigsnpr from GitHub

Source: https://github.com/privefl/bigsnpr/blob/master/README.md

Installs the bigsnpr package from its GitHub repository using the remotes package. Ensure remotes is installed first.

```r
# install.packages("remotes")
remotes::install_github("privefl/bigsnpr")
```

--------------------------------

### Example: Write and read PLINK files

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_writeBed.html

Demonstrates the process of creating a fake bigSNP object, writing it to PLINK format using snp_writeBed, and then reading it back with snp_readBed to verify data integrity. This example also shows how to attach the read data to an R session.

```R
N <- 17
M <- 911

fake <- snp_fake(N, M)
G <- fake$genotypes
G[] <- sample(as.raw(0:3), size = length(G), replace = TRUE)

# Write the object as a bed/bim/fam object
tmp <- tempfile(fileext = ".bed")
bed <- snp_writeBed(fake, tmp)

# Read this new file for the first time
rds <- snp_readBed(bed, backingfile = tempfile())
# Attach object in R session
fake2 <- snp_attach(rds)

# Same content
all.equal(fake$genotypes[], fake2$genotypes[])
#> [1] TRUE
all.equal(fake$fam, fake2$fam)
#> [1] TRUE
all.equal(fake$map, fake2$map)
#> [1] TRUE

# Two different backingfiles
fake$genotypes$backingfile
#> [1] "C:\\Users\\au639593\\AppData\\Local\\Temp\\RtmpsvW07Y\\file643850567c45.bk"
fake2$genotypes$backingfile
#> [1] "C:\\Users\\au639593\\AppData\\Local\\Temp\\RtmpsvW07Y\\file6438357c5f9b.bk"
```

--------------------------------

### Get Example Bedfile Path

Source: https://github.com/privefl/bigsnpr/blob/master/docs/articles/demo.html

Retrieves the file path for an example BED file included with the bigsnpr package. This is useful for testing and demonstration purposes.

```R
# Get the example bedfile from package bigsnpr
bedfile <- system.file("extdata", "example.bed", package = "bigsnpr")
```

--------------------------------

### Install bigsnpr package

Source: https://github.com/privefl/bigsnpr/blob/master/docs/index.html

Installation commands for the development version from GitHub or the stable version from CRAN.

```R
# install.packages("remotes")
remotes::install_github("privefl/bigsnpr")
```

```R
install.packages("bigsnpr")
```

--------------------------------

### Generate logarithmic sequences

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/seq_log.html

Examples of generating sequences with different start, end, and length parameters.

```R
seq_log(1, 1000, 4)
```

```R
seq_log(1, 100, 5)
```

--------------------------------

### Attach Example bigSNP Data

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_attachExtdata.html

Attaches an example 'bigSNP' file for use in examples and tests. This function is useful for quickly loading sample data to experiment with the package's functionalities.

```APIDOC
## Attach Example bigSNP Data

### Description
Attaches a "bigSNP" for examples and tests.

### Method
`snp_attachExtdata`

### Endpoint
N/A (R function)

### Parameters
#### Path Parameters
None

#### Query Parameters
None

#### Request Body
None

### Request Example
```R
snp_attachExtdata(bedfile = "example.bed")
```

### Arguments
- **bedfile** (character vector) - Name of one example bed file. Either `"example.bed"` (the default) or `"example-missing.bed"`.

### Value
The example "bigSNP", filebacked in the "/tmp/" directory.
```

--------------------------------

### Attach example bed file — snp_attachExtdata

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_attachExtdata.html

Use this function to attach a "bigSNP" file for examples and tests. It supports attaching either the default "example.bed" or "example-missing.bed" which contains missing data.

```R
snp_attachExtdata(bedfile = c("example.bed", "example-missing.bed"))
```

--------------------------------

### Example of GRM calculation with bed_tcrossprodSelf

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/bed_tcrossprodSelf.html

This example demonstrates how to use `bed_tcrossprodSelf` to compute a GRM from a bed file and then displays the first few rows and columns of the resulting matrix, normalized by the number of SNPs.

```R
bedfile <- system.file("extdata", "example.bed", package = "bigsnpr")
obj.bed <- bed(bedfile)

K <- bed_tcrossprodSelf(obj.bed)
K[1:4, 1:6] / ncol(obj.bed)
#>            [,1]       [,2]       [,3]       [,4]       [,5]       [,6]
#> [1,] 1.09900885 0.05531637 0.08230576 0.06440704 0.06389930 0.06632927
#> [2,] 0.05531637 1.14651667 0.08627960 0.09470694 0.12318884 0.09608808
#> [3,] 0.08230576 0.08627960 1.11142319 0.06860068 0.09951989 0.04400732
#> [4,] 0.06440704 0.09470694 0.06860068 1.06728731 0.08570288 0.09599564
    
```

--------------------------------

### Example usage of snp_ldsplit

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_ldsplit.html

Basic example showing how to load a correlation matrix for use with the splitting function.

```R
if (FALSE) {

  corr <- readRDS(url("https://www.dropbox.com/s/65u96jf7y32j2mj/spMat.rds?raw=1"))

  # adjust `THR_R2` depending on sample size used to compute corr
  # use e.g. 0.05 for small sample sizes, and 0.01 for large sample sizes
```

--------------------------------

### Create and Access a bed Object

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/sub-bed-ANY-ANY-ANY-method.html

Demonstrates how to create a 'bed' object from a file and access its elements using subsetting. Shows examples of accessing single columns, multiple columns, and checking the data type of accessed elements.

```R
bedfile <- system.file("extdata", "example-missing.bed", package = "bigsnpr")
(obj.bed <- bed(bedfile))
```

```R
obj.bed[1:5, 1]
```

```R
obj.bed[1:5, 1:2]
```

```R
typeof(obj.bed[1, 1])
```

--------------------------------

### Example Usage of bed_prodVec

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/bed_prodVec.html

Demonstrates how to use the bed_prodVec function. It first loads a bed file, creates a bed object, defines a vector y.col, and then calls bed_prodVec. The output structure is then printed.

```R
bedfile <- system.file("extdata", "example.bed", package = "bigsnpr")
obj.bed <- bed(bedfile)

y.col <- rep(1, ncol(obj.bed))
str(bed_prodVec(obj.bed, y.col))

```

--------------------------------

### Apply bed_scaleBinom to bed data

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/bed_scaleBinom.html

Examples demonstrating how to compute scaling parameters and use them within a random SVD operation.

```R
bedfile <- system.file("extdata", "example-missing.bed", package = "bigsnpr")
obj.bed <- bed(bedfile)

str(bed_scaleBinom(obj.bed))

str(bed_randomSVD(obj.bed, bed_scaleBinom))
```

--------------------------------

### Example Usage of snp_plinkQC

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_plinkQC.html

A demonstration of how to run quality control on an example dataset using the snp_plinkQC function.

```R
if (FALSE) {

bedfile <- system.file("extdata", "example.bed", package = "bigsnpr")
prefix  <- sub_bed(bedfile)
plink <- download_plink()
test <- snp_plinkQC(plink.path = plink,
                    prefix.in = prefix,
                    prefix.out = tempfile(),
                    file.type = "--bfile",  # the default (for ".bed")
                    maf = 0.05,
                    geno = 0.05,
                    mind = 0.05,
                    hwe = 1e-10,
                    autosome.only = TRUE)
test
}
```

--------------------------------

### Run PLINK IBD QC and filter samples

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_plinkIBDQC.html

This example demonstrates how to use `snp_plinkIBDQC` to compute IBD, filter related individuals based on a threshold, and save the results to a new bedfile. It also shows how to obtain the PLINK executable using `download_plink`.

```R
bedfile <- system.file("extdata", "example.bed", package = "bigsnpr")
plink <- download_plink()

bedfile <- snp_plinkIBDQC(plink, bedfile,
                          bedfile.out = tempfile(fileext = ".bed"),
                          ncores = 2)
```

--------------------------------

### Perform outlier detection with snp_pcadapt

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_pcadapt.html

Example workflow demonstrating how to attach data, perform SVD, and run the pcadapt outlier detection.

```R
test <- snp_attachExtdata()
G <- test$genotypes
obj.svd <- big_SVD(G, fun.scaling = snp_scaleBinom(), k = 10)
plot(obj.svd) # there seems to be 3 "significant" components

pcadapt <- snp_pcadapt(G, obj.svd$u[, 1:3])
snp_qq(pcadapt)
```

--------------------------------

### Match Alleles with Default Settings

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_match.html

This example demonstrates the basic usage of snp_match to align summary statistics with SNP information. It uses default settings for strand flipping and joining by position.

```R
sumstats <- data.frame(
  chr = 1,
  pos = c(86303, 86331, 162463, 752566, 755890, 758144),
  a0 = c("T", "G", "C", "A", "T", "G"),
  a1 = c("G", "A", "T", "G", "A", "A"),
  beta = c(-1.868, 0.250, -0.671, 2.112, 0.239, 1.272),
  p = c(0.860, 0.346, 0.900, 0.456, 0.776, 0.383)
)

info_snp <- data.frame(
  id = c("rs2949417", "rs115209712", "rs143399298", "rs3094315", "rs3115858"),
  chr = 1,
  pos = c(86303, 86331, 162463, 752566, 755890),
  a0 = c("T", "A", "G", "A", "T"),
  a1 = c("G", "G", "A", "G", "A")
)

snp_match(sumstats, info_snp)
```

--------------------------------

### Example Usage of snp_fastImpute

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_fastImpute.html

Demonstrates how to run the imputation, update the genotype code, and visualize the estimated error rates.

```R
if (FALSE) {

fake <- snp_attachExtdata("example-missing.bed")
G <- fake$genotypes
CHR <- fake$map$chromosome
infos <- snp_fastImpute(G, CHR)
infos[, 1:5]

# Still missing values
big_counts(G, ind.col = 1:10)
# You need to change the code of G
# To make this permanent, you need to save (modify) the file on disk
fake$genotypes$code256 <- CODE_IMPUTE_PRED
fake <- snp_save(fake)
big_counts(fake$genotypes, ind.col = 1:10)

# Plot for post-checking
## Here there is no SNP with more than 1% error (estimated)
pvals <- c(0.01, 0.005, 0.002, 0.001); colvals <- 2:5
df <- data.frame(pNA = infos[1, ], pError = infos[2, ])

# base R
plot(subset(df, pNA > 0.001), pch = 20)
idc <- lapply(seq_along(pvals), function(i) {
  curve(pvals[i] / x, from = 0, lwd = 2,
        col = colvals[i], add = TRUE)
})
legend("topright", legend = pvals, title = "p(NA & Error)",
       col = colvals, lty = 1, lwd = 2)

# ggplot2
library(ggplot2)
Reduce(function(p, i) {
  p + stat_function(fun = function(x) pvals[i] / x, color = colvals[i])
}, x = seq_along(pvals), init = ggplot(df, aes(pNA, pError))) +
  geom_point() +
  coord_cartesian(ylim = range(df$pError, na.rm = TRUE)) +
  theme_bigstatsr()
}
```

--------------------------------

### Example: PRS Calculation with GWAS, Clumping, and Thresholding

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_PRS.html

This example demonstrates a full workflow for calculating polygenic risk scores. It involves setting up data, performing PCA for covariables, conducting a GWAS on a training set, clumping significant SNPs, and then calculating PRS for a test set using various thresholds. The AUC is then plotted as a function of the number of predictors.

```R
test <- snp_attachExtdata()
G <- big_copy(test$genotypes, ind.col = 1:1000)
CHR <- test$map$chromosome[1:1000]
POS <- test$map$physical.position[1:1000]
y01 <- test$fam$affection - 1

# PCA -> covariables
obj.svd <- snp_autoSVD(G, infos.chr = CHR, infos.pos = POS)
#> Discarding 0 variant with MAC < 10 or MAF < 0.02.
#> 
#> Phase of clumping (on MAF) at r^2 > 0.2.. keep 952 variants.
#> 
#> Iteration 1:
#> Computing SVD..
#> 26 outlier variants detected..
#> 
#> Iteration 2:
#> Computing SVD..
#> 0 outlier variant detected..
#> 
#> Converged!

# train and test set
ind.train <- sort(sample(nrow(G), 400))
ind.test <- setdiff(rows_along(G), ind.train) # 117

# GWAS
gwas.train <- big_univLogReg(G, y01.train = y01[ind.train],
                                 ind.train = ind.train,
                                 covar.train = obj.svd$u[ind.train, ])
# clumping
ind.keep <- snp_clumping(G, infos.chr = CHR,
                         ind.row = ind.train,
                         S = abs(gwas.train$score))
# -log10(p-values) and thresolding
summary(lpS.keep <- -predict(gwas.train)[ind.keep])
#>     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
#> 0.000569 0.146282 0.344033 0.474463 0.640888 6.191733 
thrs <- seq(0, 4, by = 0.5)
nb.pred <- sapply(thrs, function(thr) sum(lpS.keep > thr))

# PRS
prs <- snp_PRS(G, betas.keep = gwas.train$estim[ind.keep],
               ind.test = ind.test,
               ind.keep = ind.keep,
               lpS.keep = lpS.keep,
               thr.list = thrs)

# AUC as a function of the number of predictors
aucs <- apply(prs, 2, AUC, target = y01[ind.test])
library(ggplot2)
#> Warning: package 'ggplot2' was built under R version 4.2.3
qplot(nb.pred, aucs) +

```

--------------------------------

### bed_randomSVD Function Documentation

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/bed_randomSVD.html

This section provides detailed documentation for the bed_randomSVD function, including its purpose, arguments, return values, and an example.

```APIDOC
## bed_randomSVD

### Description
Partial SVD (or PCA) of a genotype matrix stored as a PLINK (.bed) file.

### Method
`bed_randomSVD(
  obj.bed,
  fun.scaling = bed_scaleBinom,
  ind.row = rows_along(obj.bed),
  ind.col = cols_along(obj.bed),
  k = 10,
  tol = 1e-04,
  verbose = FALSE,
  ncores = 1
)`

### Parameters
#### Path Parameters
None

#### Query Parameters
None

#### Request Body
None

### Arguments
- **obj.bed** (bed-class) - Object of type [bed](bed-class.html), which is the mapping of some bed file. Use `obj.bed <- bed(bedfile)` to get this object.
- **fun.scaling** (function) - A function with parameters `X`, `ind.row` and `ind.col`, and that returns a data.frame with `$center` and `$scale` for the columns corresponding to `ind.col`, to scale each of their elements such as followed: $$\frac{X⁻\u207B - center⁻}{scale⁻}.$$ Default doesn't use any scaling. You can also provide your own `center` and `scale` by using [as_scaling_fun()](https://privefl.github.io/bigstatsr/reference/as_scaling_fun.html).
- **ind.row** (vector) - An optional vector of the row indices (individuals) that are used. If not specified, all rows are used. **Don't use negative indices.**
- **ind.col** (vector) - An optional vector of the column indices (SNPs) that are used. If not specified, all columns are used. **Don't use negative indices.**
- **k** (integer) - Number of singular vectors/values to compute. Default is `10`. **This algorithm should be used to compute only a few singular vectors/values.**
- **tol** (numeric) - Precision parameter of [svds](https://rdrr.io/pkg/RSpectra/man/svds.html). Default is `1e-4`.
- **verbose** (boolean) - Should some progress be printed? Default is `FALSE`.
- **ncores** (integer) - Number of cores used. Default doesn't use parallelism. You may use [bigstatsr::nb_cores()](https://privefl.github.io/bigstatsr/reference/reexports.html).

### Value
A named list (an S3 class "big_SVD") of:
- `d`, the singular values,
- `u`, the left singular vectors,
- `v`, the right singular vectors,
- `niter`, the number of the iteration of the algorithm,
- `nops`, number of Matrix-Vector multiplications used,
- `center`, the centering vector,
- `scale`, the scaling vector.

Note that to obtain the Principal Components, you must use [predict](https://privefl.github.io/bigstatsr/reference/predict.big_SVD.html) on the result. See examples.

### Examples
```R
bedfile <- system.file("extdata", "example.bed", package = "bigsnpr")
obj.bed <- bed(bedfile)

str(bed_randomSVD(obj.bed))
#> List of 7
#>  $ d     : num [1:10] 245.4 153.1 108.9 99.9 97.9 ...
#>  $ u     : num [1:517, 1:10] 0.0788 0.0804 0.0644 0.0777 0.0838 ...
#>  $ v     : num [1:4542, 1:10] 0.0032 -0.00161 0.02988 -0.01486 0.01259 ...
#>  $ niter : num 11
#>  $ nops  : num 182
#>  $ center: num [1:4542] 0.685 0.412 0.474 0.369 0.913 ...
#>  $ scale : num [1:4542] 0.671 0.572 0.601 0.549 0.704 ...
#>  - attr(*, "class")= chr "big_SVD"
```
```

--------------------------------

### Install bigsnpr from CRAN

Source: https://github.com/privefl/bigsnpr/blob/master/README.md

Installs the stable version of the bigsnpr package directly from the Comprehensive R Archive Network (CRAN).

```r
install.packages("bigsnpr")
```

--------------------------------

### Example usage of snp_MAX3

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_MAX3.html

Demonstrates generating fake genotype data and calculating MAX3 statistics, including the Armitage trend test variant.

```R
set.seed(1)

# constructing a fake genotype big.matrix
N <- 50; M <- 1200
fake <- snp_fake(N, M)
G <- fake$genotypes
G[] <- sample(as.raw(0:3), size = length(G), replace = TRUE)
G[1:8, 1:10]
#>      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
#> [1,]    0    0   NA    0    2    1    0   NA   NA     2
#> [2,]   NA    0   NA    2   NA    2    2    2    2     1
#> [3,]    2    2    2    2    1    2    1    0    2     0
#> [4,]   NA    1    0   NA    0    2    2    0    0     2
#> [5,]    0    0   NA   NA    2   NA    2    2   NA    NA
#> [6,]    1    0   NA    2    2   NA    1    2    2     0
#> [7,]    2   NA    1    1    1   NA    1    2    1    NA
#> [8,]    2    2    1    0    1   NA   NA    2    0     1

# Specify case/control phenotypes
fake$fam$affection <- rep(1:2, each = N / 2)

# Get MAX3 statistics
y01 <- fake$fam$affection - 1
str(test <- snp_MAX3(fake$genotypes, y01.train = y01))
#> Classes 'mhtest' and 'data.frame':	1200 obs. of  1 variable:
#>  $ score: num  1.7124 0.0178 2.7611 0.1309 0.4396 ...
#>  - attr(*, "transfo")=function (x)  
#>  - attr(*, "predict")=function (xtr)  
# p-values are not well calibrated
snp_qq(test)

# genomic control is not of much help
snp_qq(snp_gc(test))


# Armitage trend test (well calibrated because only one test)
test2 <- snp_MAX3(fake$genotypes, y01.train = y01, val = 0.5)
snp_qq(test2)
```

--------------------------------

### Generate Manhattan and QQ plots

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_manhattan.html

Example workflow for performing a GWAS, applying genomic control, and visualizing results with QQ and Manhattan plots.

```R
set.seed(9)

test <- snp_attachExtdata()
G <- test$genotypes
y <- rnorm(nrow(G))

gwas <- big_univLinReg(G, y)

snp_qq(gwas)

gwas_gc <- snp_gc(gwas) # this modifies `attr(gwas_gc, "transfo")`
snp_qq(gwas_gc)


# The next plot should be prettier with a real dataset
snp_manhattan(gwas_gc,
              infos.chr = test$map$chromosome,
              infos.pos = test$map$physical.pos) +
  ggplot2::geom_hline(yintercept = -log10(5e-8), linetype = 2, color = "red")


p <- snp_qq(gwas_gc) +
  ggplot2::aes(text = asPlotlyText(test$map)) +
  ggplot2::labs(subtitle = NULL, x = "Expected -log10(p)", y = "Observed -log10(p)")
if (FALSE) plotly::ggplotly(p, tooltip = "text")
```

--------------------------------

### Visualize LD-split trade-offs using ggplot2

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_ldsplit.html

These examples use `ggplot2` to visualize the trade-offs between different parameters (number of blocks, percentage kept, sum of squared sizes) and the cost associated with LD splitting. The plots help in selecting an optimal `max_size`.

```R
library(ggplot2)
# trade-off cost / number of blocks
qplot(n_block, cost, color = factor(max_size, SEQ), data = res) + 
  theme_bw(14) + 
  scale_y_log10() + 
  theme(legend.position = "top") + 
  labs(x = "Number of blocks", color = "Maximum block size",
       y = "Sum of squared correlations outside blocks")
```

```R
# trade-off cost / number of non-zero values
qplot(perc_kept, cost, color = factor(max_size, SEQ), data = res) + 
  theme_bw(14) + 
  # scale_y_log10() + 
  theme(legend.position = "top") + 
  labs(x = "Percentage of non-zero values kept", color = "Maximum block size",
       y = "Sum of squared correlations outside blocks")
```

```R
# trade-off cost / sum of squared sizes
qplot(cost2, cost, color = factor(max_size, SEQ), data = res) + 
  theme_bw(14) + 
  scale_y_log10() + 
  geom_vline(xintercept = 0)+
  theme(legend.position = "top") + 
  labs(x = "Sum of squared blocks", color = "Maximum block size",
       y = "Sum of squared correlations outside blocks")
```

--------------------------------

### Match Alleles without Strand Flip

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_match.html

This example shows how to use snp_match while disabling the strand flip functionality. This is useful when you are certain about the strand orientation or want to handle it manually.

```R
sumstats <- data.frame(
  chr = 1,
  pos = c(86303, 86331, 162463, 752566, 755890, 758144),
  a0 = c("T", "G", "C", "A", "T", "G"),
  a1 = c("G", "A", "T", "G", "A", "A"),
  beta = c(-1.868, 0.250, -0.671, 2.112, 0.239, 1.272),
  p = c(0.860, 0.346, 0.900, 0.456, 0.776, 0.383)
)

info_snp <- data.frame(
  id = c("rs2949417", "rs115209712", "rs143399298", "rs3094315", "rs3115858"),
  chr = 1,
  pos = c(86303, 86331, 162463, 752566, 755890),
  a0 = c("T", "A", "G", "A", "T"),
  a1 = c("G", "G", "A", "G", "A")
)

snp_match(sumstats, info_snp, strand_flip = FALSE)
```

--------------------------------

### Get sample information from external files

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_getSampleInfos.html

Use `snp_getSampleInfos` to retrieve information of individuals by matching from an external file. This function can take a data.frame or a character vector of file names as input. Ensure the external files contain columns for family and sample IDs.

```R
test <- snp_attachExtdata()
table(test$fam$family.ID)
```

```R
files <- system.file("extdata", paste0("cluster", 1:3), package = "bigsnpr")
bigreadr::fread2(files[1])
```

--------------------------------

### Initialize Environment and Data

Source: https://github.com/privefl/bigsnpr/blob/master/docs/articles/bedpca.html

Load necessary libraries and download the 1000 Genomes subset data.

```R
library(ggplot2)
library(bigsnpr)

bedfile <- download_1000G("tmp-data")
```

--------------------------------

### Download and Unzip Tutorial Data

Source: https://github.com/privefl/bigsnpr/blob/master/docs/articles/LDpred2.html

Downloads a ZIP archive of tutorial data and extracts it into a local directory.

```R
# install.packages("runonce")
zip <- runonce::download_file(
  "https://github.com/privefl/bigsnpr/raw/master/data-raw/public-data3.zip",
  dir = "tmp-data")
unzip(zip)
```

--------------------------------

### Constructor: bed(bedfile)

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/bed-class.html

Initializes a new bed object by providing the path to a .bed file. The corresponding .bim and .fam files must exist in the same directory.

```APIDOC
## Constructor: bed(bedfile)

### Description
Creates a reference class object that points to a mapped version of a PLINK bed file.

### Parameters
#### Arguments
- **bedfile** (string) - Required - Path to the .bed file. Requires corresponding .bim and .fam files in the same directory.

### Request Example
bedfile <- system.file("extdata", "example-missing.bed", package = "bigsnpr")
obj.bed <- bed(bedfile)
```

--------------------------------

### Download and Prepare Genotype Data

Source: https://github.com/privefl/bigsnpr/blob/master/docs/articles/ancestry.html

Downloads BED, FAM, and BIM files from a URL and unzips the BIM file. Ensure the directory exists before running.

```R
prefix <- "https://sharehost.hms.harvard.edu/genetics/reich_lab/sgdp/variant_set/cteam_extended.v4.maf0.1perc"
bedfile <- runonce::download_file(paste0(prefix, ".bed"), dir = DIR)
famfile <- runonce::download_file(paste0(prefix, ".fam"), dir = DIR)
bim_zip <- runonce::download_file(paste0(prefix, ".bim.zip"), dir = DIR)
unzip(bim_zip, exdir = DIR, overwrite = FALSE)
```

--------------------------------

### Initialize Polygenic Risk Scores (PRS) Workflow

Source: https://context7.com/privefl/bigsnpr/llms.txt

Sets up training and test data splits for PRS calculation using bigsnpr.

```r
library(bigsnpr)

# Load example data
snp_data <- snp_attachExtdata()
G <- big_copy(snp_data$genotypes, ind.col = 1:1000)
CHR <- snp_data$map$chromosome[1:1000]
POS <- snp_data$map$physical.pos[1:1000]
y <- snp_data$fam$affection - 1

# Split into training and test sets
ind.train <- sample(nrow(G), 400)
ind.test <- setdiff(1:nrow(G), ind.train)
```

--------------------------------

### GET /download_genetic_map

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_asGeneticPos2.html

Downloads a specified genetic map file to a local directory.

```APIDOC
## GET /download_genetic_map

### Description
Downloads genetic map data for specified genome builds (hg19 or hg38) to a local directory.

### Method
GET

### Parameters
#### Query Parameters
- **type** (string) - Required - The genetic map type: "hg19_OMNI", "hg19_hapmap", or "hg38_price".
- **dir** (string) - Required - Directory path where files will be downloaded and decompressed.
- **ncores** (integer) - Optional - Number of cores to use for parallel processing. Default is 1.
```

--------------------------------

### Create a bed object

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/bed-class.html

Initializes a bed object by providing the path to a .bed file. Requires corresponding .bim and .fam files in the same directory. This object stores a pointer to the mapped bed file.

```R
bedfile <- system.file("extdata", "example-missing.bed", package = "bigsnpr")
(obj.bed <- bed(bedfile))
```

--------------------------------

### Example usage of snp_cor

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_cor.html

Demonstrates calculating a correlation matrix for a subset of SNPs from a genotype object.

```R
test <- snp_attachExtdata()
G <- test$genotypes

corr <- snp_cor(G, ind.col = 1:1000)
corr[1:10, 1:10]
```

--------------------------------

### Generate Manhattan Plot

Source: https://github.com/privefl/bigsnpr/blob/master/docs/articles/demo.html

Generates a Manhattan plot for a GWAS object. This example assumes data for a single chromosome.

```r
# Manhattan plot, not quite sexy because there are only 1 chromosome here
snp_manhattan(obj.gwas.gc, infos.chr = CHR, infos.pos = POS)
```

--------------------------------

### Customize and visualize GWAS plots

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_qq.html

Examples of creating Manhattan plots and interactive Q-Q plots using ggplot2 and plotly extensions.

```R
# The next plot should be prettier with a real dataset
snp_manhattan(gwas_gc,
              infos.chr = test$map$chromosome,
              infos.pos = test$map$physical.pos) +
  ggplot2::geom_hline(yintercept = -log10(5e-8), linetype = 2, color = "red")


p <- snp_qq(gwas_gc) +
  ggplot2::aes(text = asPlotlyText(test$map)) +
  ggplot2::labs(subtitle = NULL, x = "Expected -log10(p)", y = "Observed -log10(p)")
if (FALSE) plotly::ggplotly(p, tooltip = "text")
```

--------------------------------

### Solve Quadratic Programming for Projected Individuals

Source: https://github.com/privefl/bigsnpr/blob/master/docs/articles/ancestry.html

Uses quadprog::solve.QP to compute projections for each individual, aggregating results by group factor.

```R
bvec <- c(1, rep(0, ncol(X)))

# solve a QP for each projected individual
all_res <- apply(all_proj, 1, function(y) {
  quadprog::solve.QP(
    Dmat = cp_X_pd$mat,
    dvec = crossprod(y, X),
    Amat = Amat,
    bvec = bvec,
    meq  = 1
  )$sol %>% 
    tapply(grp_fct, sum) %>% 
    round(7)
})
```

--------------------------------

### Calculate allele frequencies for a bed file

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/bed_MAF.html

Example usage of bed_MAF to calculate allele frequencies for a subset of SNPs from a bed file.

```R
bedfile <- system.file("extdata", "example-missing.bed", package = "bigsnpr")
obj.bed <- bed(bedfile)

bed_MAF(obj.bed, ind.col = 1:5)
```

--------------------------------

### Prepare for Quadratic Programming

Source: https://github.com/privefl/bigsnpr/blob/master/docs/articles/ancestry.html

Calculates the cross-product of the projected reference allele frequencies and finds the nearest positive-definite matrix. Sets up the Amat matrix for quadratic programming.

```R
cp_X_pd <- Matrix::nearPD(crossprod(X), base.matrix = TRUE)
Amat <- cbind(1, diag(ncol(X)))
```

--------------------------------

### Run LDpred2-grid and Evaluate

Source: https://github.com/privefl/bigsnpr/blob/master/docs/articles/LDpred2.html

Execute the LDpred2-grid algorithm and calculate Z-scores for each model configuration.

```R
set.seed(1)  # to get the same result every time
# takes less than 2 min with 4 cores
beta_grid <- snp_ldpred2_grid(corr, df_beta, params, ncores = NCORES)

pred_grid <- big_prodMat(G, beta_grid, ind.col = df_beta[["_NUM_ID_"]])
params$score <- apply(pred_grid[ind.val, ], 2, function(x) {
  if (all(is.na(x))) return(NA)
  summary(lm(y[ind.val] ~ x))$coef["x", 3]
  # summary(glm(y[ind.val] ~ x, family = "binomial"))$coef["x", 3]  # for a binary phenotype
})
```

--------------------------------

### Plot Manhattan with Genomic Control

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_gc.html

Generates a Manhattan plot from GWAS results after applying Genomic Control. This example adds a horizontal line for the significance threshold.

```R
snp_manhattan(gwas_gc,
                  infos.chr = test$map$chromosome,
                  infos.pos = test$map$physical.pos) + 
  ggplot2::geom_hline(yintercept = -log10(5e-8), linetype = 2, color = "red")
```

--------------------------------

### Read UK Biobank BGEN files

Source: https://context7.com/privefl/bigsnpr/llms.txt

Use snp_readBGEN() to read UK Biobank BGEN files. Ensure you have the required packages installed for BGEN file handling.

```r
obj.bigSNP <- snp_readBGEN("path/to/your/data.bgen")
```

--------------------------------

### Perform partial SVD on a bed file

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/bed_randomSVD.html

Initializes a bed object from a file and computes the partial SVD using default parameters.

```R
bedfile <- system.file("extdata", "example.bed", package = "bigsnpr")
obj.bed <- bed(bedfile)

str(bed_randomSVD(obj.bed))
```

--------------------------------

### Visualize chain convergence

Source: https://github.com/privefl/bigsnpr/blob/master/docs/articles/LDpred2.html

Plots the path of p and h2 estimates for a chain to verify convergence.

```R
library(ggplot2)
auto <- multi_auto[[1]]  # first chain
plot_grid(
  qplot(y = auto$path_p_est) + 
    theme_bigstatsr() + 
    geom_hline(yintercept = auto$p_est, col = "blue") +
    scale_y_log10() +
    labs(y = "p"),
  qplot(y = auto$path_h2_est) + 
    theme_bigstatsr() + 
    geom_hline(yintercept = auto$h2_est, col = "blue") +
    labs(y = "h2"),
  ncol = 1, align = "hv"
)
```

--------------------------------

### Calculate Fst between SNP sets

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_fst.html

Calculates the Fst statistic between specified sets of SNPs. Use `overall = TRUE` to get a single Fst value for the entire set.

```R
snp_fst(list_df_af[c(1, 2)], overall = TRUE)
```

```R
snp_fst(list_df_af[c(1, 3)], overall = TRUE)
```

```R
snp_fst(list_df_af[c(3, 2)], overall = TRUE)
```

--------------------------------

### Select and visualize blocks from LD-split results

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_ldsplit.html

This snippet shows how to select an optimal split based on a cost-performance metric and then transform the sparse correlation matrix into a triplet format for visualization. It requires the `dplyr` package.

```R
library(dplyr)
all_ind <- res %>%
  arrange(cost2 * sqrt(5 + cost)) %>%
  print() %>%
  slice(1) %>%
  pull(all_last)
```

```R
## Transform sparse representation into (i,j,x) triplets
corrT <- as(corr, "dgTMatrix")
upper <- (corrT@i <= corrT@j & corrT@x^2 >= THR_R2)
df <- data.frame(
  i = corrT@i[upper] + 1L,
  j = corrT@j[upper] + 1L,
  r2 = corrT@x[upper]^2
)
df$y <- (df$j - df$i) / 2
```

```R
ggplot(df) + 
  geom_point(aes(i + y, y, alpha = r2)) + 
  theme_minimal() + 
  theme(axis.text.y = element_blank(), axis.ticks.y = element_blank(),
        strip.background = element_blank(), strip.text.x = element_blank()) + 
  scale_alpha_continuous(range = 0:1) + 
  scale_x_continuous(expand = c(0.02, 0.02), minor_breaks = NULL,
                     breaks = head(all_ind[[1]], -1) + 0.5) + 
  facet_wrap(~ cut(i + y, 4), scales = "free", ncol = 1) + 
  labs(x = "Position", y = NULL)
```

--------------------------------

### Example of snp_ldsc2 usage

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_ldsc.html

Demonstrates how to use the `snp_ldsc2` function with a correlation matrix and GWAS beta estimates. The first call estimates both intercept and heritability, while the second call fixes the intercept to 1.

```R
bigsnp <- snp_attachExtdata()
    G <- bigsnp$genotypes
    y <- bigsnp$fam$affection - 1
    corr <- snp_cor(G, ind.col = 1:1000)
    
    gwas <- big_univLogReg(G, y, ind.col = 1:1000)
    df_beta <- data.frame(beta = gwas$estim, beta_se = gwas$std.err,
                          n_eff = 4 / (1 / sum(y == 0) + 1 / sum(y == 1)))
    
    snp_ldsc2(corr, df_beta)
    #>       int        h2 
    #> 1.0000000 0.2335429 
    snp_ldsc2(corr, df_beta, blocks = 20, intercept = NULL)
    #>       int    int_se        h2     h2_se 
    #> 0.4986445 0.2526338 0.6195226 0.1818980
```

--------------------------------

### Project PCA and Visualize

Source: https://github.com/privefl/bigsnpr/blob/master/docs/articles/bedpca.html

Projects individuals onto existing PCA components and visualizes the results. Ensure necessary objects like obj.svd2, obj.bed, and PCs are defined. The projection might be 'useless' if too few individuals are projected.

```R
proj <- bed_projectSelfPCA(obj.svd2, obj.bed, 
                               ind.row = rows_along(obj.bed)[-ind.row],
                               ncores = 1) # useless -> too few individuals
    PCs[-ind.row, ] <- proj$OADP_proj

    plot(PCs[ind.row, 7:8], pch = 20, xlab = "PC7", ylab = "PC8")
    points(PCs[-ind.row, 7:8], pch = 20, col = "blue")
```

--------------------------------

### Create Interactive QQ Plot with Genomic Control

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_gc.html

Creates an interactive QQ plot using `plotly::ggplotly` after applying Genomic Control. This example customizes plot labels and adds text for tooltips.

```R
p <- snp_qq(gwas_gc) +
  ggplot2::aes(text = asPlotlyText(test$map)) +
  ggplot2::labs(subtitle = NULL, x = "Expected -log10(p)", y = "Observed -log10(p)")
if (FALSE) plotly::ggplotly(p, tooltip = "text")
```

--------------------------------

### Get Posterior Beta Samples with LDpred2

Source: https://github.com/privefl/bigsnpr/blob/master/docs/articles/prs_uncertainty.html

Obtain MCMC samples from the posterior of causal effects by supplying tuned hyperparameters to snp_ldpred2_grid with return_sampling_betas = TRUE. Ensure best_param is a data frame with exactly one row.

```R
best_param <- data.frame(p = 0.01, h2 = 0.2, sparse = FALSE)

posterior_beta_samples <- snp_ldpred2_grid(
  corr, df_beta, best_param,
  return_sampling_betas = TRUE, num_iter = 500)
dim(posterior_beta_samples)
```

--------------------------------

### Remove samples based on IBD data

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_plinkIBDQC.html

This example shows how to use `snp_plinkRmSamples` to remove individuals from a bedfile based on IBD data. It filters the `df_rel` data frame to include only individuals with `PI_HAT > 0.2` and passes this subset to the function.

```R
snp_plinkRmSamples(plink, bedfile,
                       bedfile.out = tempfile(fileext = ".bed"),
                       df.or.files = subset(df_rel, PI_HAT > 0.2))
```

--------------------------------

### Retrieve and summarize sample information

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_getSampleInfos.html

Extracts sample information from files and displays a frequency table of the populations.

```R
infos <- snp_getSampleInfos(test, files, header = FALSE)
table(infos[[1]])
```

--------------------------------

### Visualize IBD data

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_plinkIBDQC.html

These examples demonstrate how to visualize the IBD data obtained from `snp_plinkIBDQC` using `ggplot2`. The first plot visualizes the relationship between Z0 and Z1, colored by RT, while the second plots PI_HAT with a horizontal line indicating a threshold.

```R
library(ggplot2)
qplot(Z0, Z1, data = df_rel, col = RT)
```

```R
qplot(y = PI_HAT, data = df_rel) +
  geom_hline(yintercept = 0.2, color = "blue", linetype = 2)
```

--------------------------------

### Download PLINK 2.0 Executable

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/download_plink.html

Downloads the PLINK 2.0 executable. Allows selection of AVX2, ARM, or AMD optimized versions, along with directory, overwrite, and verbosity options.

```R
download_plink2(
  dir = tempdir(),
  AVX2 = TRUE,
  ARM = FALSE,
  AMD = FALSE,
  overwrite = FALSE,
  verbose = TRUE
)
```

--------------------------------

### Get IBD data without blind QC

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/snp_plinkIBDQC.html

This snippet shows how to run `snp_plinkIBDQC` with `do.blind.QC = FALSE` to obtain a data frame of IBD information without performing automatic sample filtering. This allows for manual inspection and visualization of relatedness.

```R
df_rel <- snp_plinkIBDQC(plink, bedfile, do.blind.QC = FALSE, ncores = 2)
str(df_rel)
```

--------------------------------

### Retrieve dimensions and length for bed and bed_light objects

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/bed-methods.html

These methods provide the dimensions and length of bed and bed_light objects. They are automatically used by nrow() and ncol().

```R
# S4 method for bed
dim(x)

# S4 method for bed
length(x)

# S4 method for bed_light
dim(x)

# S4 method for bed_light
length(x)
```

--------------------------------

### Download PLINK 1.9 Executable

Source: https://github.com/privefl/bigsnpr/blob/master/docs/reference/download_plink.html

Downloads the PLINK 1.9 executable. Specify the directory for download and whether to overwrite existing files or show verbose output.

```R
download_plink(dir = tempdir(), overwrite = FALSE, verbose = TRUE)
```