1
$\begingroup$

I have a (SNP) distance matrix that I would like to convert to coordinates (ideally just one linear scale) in R. I am planning to use the transformed coordinates in principle component analysis along with other variables and would also like to plot the coordinates against e.g. time.

The S/O post here gives an example that seems to do what I want, with multi-dimensional scaling. It is also the approach taken in the bioconductor package SNPrelate (see here) However when I test this out on some dummy data, the distance between the coordinate points does not seem to correspond at all to the original distance defined by the matrix and I'm not sure what I am doing wrong / missing.

Here is a MWE:

# Load libraries
library(data.table)
library(MASS)
library(ggplot2)
library(plotly)

# Function to create dummy matrix of n size with y dummy values:
# Derived from: https://davetang.org/muse/2014/01/22/making-symmetric-matrices-in-r/

CreateMatrix <- function(size, minnumber, maxnumber) {

  # Calculate the number of combinations:
  ncombos = choose(size, 2)

  # Create the sequence of dummy values:
  values = sample(c(minnumber:maxnumber), replace = T, size = ncombos)

  # Create an empty matrix of 0s with rows and columns labelled with capital letters:
  mymatrix = matrix(rep(0, size*size), nrow = size, dimnames = list(LETTERS[1:size]))

  # Fill the the matrix with sequence of dummy values:
  mymatrix[lower.tri(mymatrix)] = values
  mymatrix = t(mymatrix)
  mymatrix[lower.tri(mymatrix)] = values

  # Add row names:
  rownames(mymatrix) = c(LETTERS[1:size])

  # Return the completed matrix:
  return(mymatrix)

}

############################################################################
# CREATE THE DUMMY MATRIX:

# Create a dummy matrix for:
# - SNP distances between 10 isolates
# - where the minimum SNP distance is 0
# - and the maximum SNP distance is 20

snpmat <- CreateMatrix(size = 10, minnumber = 0, maxnumber = 20)

############################################################################

Then converting to a single dimension scale (I also tried with 2 dimensions but it gives the same problem) using cmdscale from the stats package:

# GET X (AND Y IF DESIRED) COORDINATES FOR MATRIX VALUES:

# Use multi-dimensional scaling to get x and y coordinates from the SNP distances:
mds <- cmdscale(snpmat, k = 1, add = F) # k = number of dimensions, default is 2 (x and y)

# Convert the results into a data.table with an ID column for merging on:
mds <- data.table(id = dimnames(mds)[[1]], mds)

# Give the coordinates nice names:
setnames(mds, old = c("V1"), new = c("snp_x"))

... and plotting the results to view:

# Create a scatter plot of the results and check if sensible:
snpdist <- ggplot(data = mds, aes(snp_x, y = 1)) +
  geom_point(aes(colour = id))

# View the SNP distance matrix:
snpmat

# View the plot:
snpdist

# Convert to interactive for sense-check:
ggplotly()

Isolates that should be 0 SNPs apart according to the matrix are being represented at various distances from each other as a coordinate, while I would have expected them to overlap completely on the above plot. In general the coordinates bear no relation to the original distance matrix (that I can see).

$\endgroup$
  • $\begingroup$ I am now considering other approaches, such as expressing each pairwise distance as a proportion of the maximum pairwise distance in my data set - but would still like to understand what multi-dimensional scaling is doing in the above example (and why it is not doing what I thought it should do). $\endgroup$ – Amy M Jun 26 '18 at 17:48

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.