library(data.table)
library(fasttime)
library(ggplot2)
library(PerformanceAnalytics)
library(xts)

#Read data
spy = fread("/path/to/unzipped/tsv")

#Convert dates to POSIX format
spy = spy[,date := as.POSIXct(strptime(date, format="%Y-%m-%d %H:%M:%S"), tz="UTC")]

#Convert to time series object
spy_xts = as.xts.data.table(spy)

#Convert timezone to EST
indexTZ(spy_xts) <- "America/New_York"

#Select data points only at 9:30, 10:09, and 10:20 AM
spy_xts =rbind(spy_xts['T09:30/T09:31'], spy_xts['T10:09/T10:11'], spy_xts['T10:20/T10:21'])

#Calculate interval returns based on close price
spy_xts$Returns = Return.calculate(spy_xts$close, method="discrete")
spy_xts$Log_Returns = Return.calculate(spy_xts$close, method="log")

#Make a matrix of returns at 10:09 and 10:20
rets = cbind(as.vector(spy_xts['T10:09/T10:11']$Log_Returns), as.vector(spy_xts['T10:20/T10:21']$Log_Returns))

#Scatter plot of returns at 10:09 vs 10:20
plot(rets)

#Pearson and Spearman correlations
cor(rets)
cor(rets, method="spearman")