Standard Deviation in R: sd(), na.rm, and By-Group Calculation

sd(x) computes the sample standard deviation in R — but if x contains even one NA, sd() silently returns NA for the whole vector unless you add na.rm = TRUE. There is no warning and no error; the function just returns a missing value, and it is easy to mistake that for “no variation in the data” rather than “some values were missing.” This page covers that gotcha, what R’s sd() actually calculates (the sample, not the population, formula), and the two standard ways to get standard deviation by group: base R’s aggregate() and dplyr::group_by() with summarise().

The basic call: sd(x)

sd() is a base R function — no package needed. Its signature is:

sd(x, na.rm = FALSE)

Called on a plain numeric vector with no missing values, it just works:

scores <- c(78, 85, 90, 88, 76, 95, 89) sd(scores) #> [1] 6.998299

R’s sd() always computes the sample standard deviation — it divides the sum of squared deviations by n - 1 (Bessel’s correction), not by n. There is no built-in argument or companion function to get the population standard deviation (dividing by n) directly; if you specifically need that, compute it by hand:

pop_sd <- function(x, na.rm = FALSE) {   if (na.rm) x <- x[!is.na(x)]   sqrt(sum((x - mean(x))^2) / length(x)) }

For almost all research-data use cases — a sample of participants, trials, or observations standing in for a larger population — the sample formula sd() already gives you is the correct one to report; population SD only applies when your data genuinely is the entire population you care about, which is rare outside a full-census dataset.

The na.rm trap

This is the part that actually catches people. Add one missing value to the vector above and the call changes behavior completely:

scores_na <- c(78, 85, 90, NA, 76, 95, 89) sd(scores_na) #> [1] NA

No warning, no error — just NA. This is consistent with how most of R’s summary statistics functions behave (mean(), sum(), var(), median() all default to na.rm = FALSE too), but it is a common source of silently wrong output in a pipeline: a script that runs sd() inside a loop or a summarise() call over many columns or groups can return a column of NAs with no indication anything went wrong, especially if you are not eyeballing every value. The fix is one argument:

sd(scores_na, na.rm = TRUE) #> [1] 7.359801

With na.rm = TRUE, the missing value is dropped before the mean and the squared deviations are computed — note the result (7.36) is not the same number you’d get by treating the missing value as anything other than absent; it is the SD of the six remaining non-missing values. Two things worth checking whenever you use na.rm = TRUE in a real analysis:

• Report how many values were dropped. sum(is.na(scores_na)) tells you the count. If missingness is substantial or non-random (e.g., concentrated in one condition or one site), silently dropping it and reporting only the SD of what’s left can misrepresent the data — that’s a methods-reporting issue, not just a coding one.
• na.rm does not distinguish why a value is missing. R only sees NA; it has no concept of missing-at-random versus missing-not-at-random. Statistical treatment of missingness (imputation, sensitivity analysis) is a separate decision from the na.rm argument, which is purely “include this row in the arithmetic or don’t.”

Standard deviation by group

Most real analyses need SD broken out by a grouping variable — by condition, by site, by year — not one number for an entire dataset. R has two standard routes to this: base R’s aggregate(), and the tidyverse’s dplyr::group_by() plus summarise(). Both give the same numbers; they differ in syntax and in what else they make convenient.

Example data

df <- data.frame(   group = c("A", "A", "A", "B", "B", "B", "B"),   value = c(12, 15, 11, 22, NA, 19, 25) )

Base R: aggregate()

aggregate() uses a formula interface: value ~ group reads as “compute this for value, split by group.” The function you want to apply goes in the FUN argument, and any extra arguments that function needs — like na.rm — are passed straight through after it:

aggregate(value ~ group, data = df, FUN = sd, na.rm = TRUE) #>   group    value #> 1     A 2.081666 #> 2     B 3.214550

One easy-to-miss detail: if you omit na.rm = TRUE here, aggregate()‘s formula interface drops entire rows with an NA in any variable in the formula before the calculation even runs (its default na.action is na.omit) — so group B’s row with the missing value disappears from the input entirely, rather than sd() returning NA for that group. The output looks clean either way, which is exactly why it’s worth checking your group sizes (table(df$group)) rather than assuming the row count matches what you expect.

Tidyverse: dplyr::group_by() + summarise()

library(dplyr)  df %>%   group_by(group) %>%   summarise(sd_value = sd(value, na.rm = TRUE), n = n()) #> # A tibble: 2 x 3 #>   group sd_value     n #>   <chr>    <dbl> <int> #> 1 A         2.08     3 #> 2 B         3.21     4 

Here na.rm = TRUE is passed directly to sd() the same way it would be outside a pipe, and unlike aggregate()‘s formula interface, group_by()/summarise() does not drop the whole row for group B — it keeps all four rows in the group and n still reports 4, while sd() itself only excludes the one missing value from its own calculation. Adding n = n() alongside the SD is good practice generally: a standard deviation reported without its group size is hard to interpret, and it makes a shrunken group (from missingness or an unbalanced design) visible in the same table rather than hidden.

summarise() also composes cleanly with multiple statistics and multiple grouping variables in one call:

df %>%   group_by(group) %>%   summarise(     mean_value = mean(value, na.rm = TRUE),     sd_value = sd(value, na.rm = TRUE),     n = n()   )

aggregate() vs. dplyr: which to use

• Base R (aggregate()) needs no package install, which matters for a script meant to run anywhere with only base R, or a course context where installing packages is friction. Its formula interface is compact for one or two summary statistics, but its silent row-dropping behavior on missing data (above) is worth knowing before you trust a clean-looking result.
• dplyr is more explicit about what happens to missing values (it only removes what you tell sd() to remove, not the whole row), reads more clearly once you’re computing several statistics at once, and is the more common convention in current published R analysis code and teaching material. It requires library(dplyr) (part of the tidyverse package collection).

Neither is “more correct” — they compute identical numbers when used carefully. The practical difference that matters most for reproducibility is the missing-data handling difference above: know which one your script is doing, and report your final group sizes rather than assuming they match your raw data’s row counts.

Reporting standard deviation correctly

A standard deviation reported without stating n (or without stating that it’s the sample SD, which is the R default) is incomplete for reproducibility purposes — a reader checking your work needs both numbers to reconstruct the calculation, and needs to know how missing data was handled to reproduce your exact figure. This is the same reproducibility standard that applies to any descriptive statistic reported in a manuscript or a data management plan appendix: state the formula convention used (sample vs. population), the software and function, and how missing values were treated, not just the resulting number.

Related reading

• Standard Deviation: Formula, Worked Example, and How to Interpret It — the underlying statistic itself, independent of any software.
• Running ANOVA in R: aov(), car::Anova(), and Type III Sums of Squares — another common R gotcha (default sums-of-squares type) in the same statistical-software track.
• Logistic Regression in R: glm(), Odds Ratios, and Diagnostics

Frequently asked questions

Why does R’s sd() give a different number than a calculator or textbook formula?

Usually because of the denominator: R’s sd() always divides by n - 1 (the sample formula, Bessel’s correction), while some basic calculators or intro textbook examples use n (the population formula). For a sample size of, say, 10, the two will differ by roughly 5%. Check which formula the comparison source is using before assuming sd() is wrong.

Does na.rm = TRUE affect the mean used inside the SD calculation, or just which values are counted?

Both, and consistently: with na.rm = TRUE, sd() first drops the NA values, then computes the mean of what remains, then computes the sum of squared deviations from that same mean, divided by (remaining n) − 1. It is not mixing a mean computed on the full vector with an SD computed on the reduced one.

Can I get standard deviation for multiple numeric columns at once, by group?

Yes, with dplyr::across() inside summarise(): df %>% group_by(group) %>% summarise(across(where(is.numeric), (x) sd(x, na.rm = TRUE))) applies sd() to every numeric column in one call. Base R’s equivalent is aggregate(cbind(col1, col2) ~ group, data = df, FUN = sd, na.rm = TRUE).

Is there a function for standard error of the mean in base R?

No built-in se() function exists in base R. It’s computed directly as sd(x, na.rm = TRUE) / sqrt(sum(!is.na(x))) — the sample SD divided by the square root of the number of non-missing observations.