Having trouble conceptualising how to implement a linear interpolation in R

Question

I am fairly new to coding and am having trouble even conceptualising how I will implement this problem in R. Sorry if I am missing something totally obvious.

I am trying to find the pressure correction factor (y) for each minute at a particular meteorological station. This is the formula:

x is the ratio of the elevation of the location over some constant

y is the pressure correction factor

m is the optical air mass at sea level

m depends on the time of day, and my time series has a column for the value of m at each minute of the series.

when x = 1, y = 1.
when x = 0.75, y= 1.25 -.012*m  + .00037*m^2
when x= 0.5, y = 1.68 -.03*m + .0009*m^2

Between values of x, y is interpolated linearly.

The value of x I am interested in is 0.9804

This is as far as I have gotten in terms of implementation (it is not far):

x <- c(0.5 , .75 , 1)
pc1 = 1.248274 - 0.011997*DF$AirMass0 + 0.000370*DF$AirMass0^2
pc2 = 1.68219 - 0.03059*DF$AirMass0 + 0.000890*DF$AirMass0^2

pc <- c(0 , pc1[1], pc2[1] , 1)
approx(x, pc , n = 1000)

I know this is just calling my first value of airmass in my time series. Also, that even with n=1000, I will not be able to get an exact value of x=9804.

I realise it is a long shot that anyone will even read all of this, but I would really appreciate any help.

Thanks in advance :)

EDIT: Here's the head of my df

> dput(head(DF2))
structure(list(DF.Date = structure(c(1357296180, 1358021760, 
1357991280, 1357629540, 1357287540, 1358158260), class = c("POSIXct", 
"POSIXt"), tzone = "UTC"), DF.Irradiance = c(1055.64, 0, 975.22, 
0.25, 953.33, 1017.62), DF.AirMass0 = c(0.0279842300461168, 0.0195227712650674, 
0.0194685363977006, 0.0384509898910604, -0.142855127416022, 0.0322758934562728
)), row.names = c(NA, 6L), class = "data.frame")


> StationHeightCorrection <- exp(-167/8435.52)
> StationHeightCorrection
[1] 0.9803974

Could you also provide `DF`, at least part of it? You could try `dput(head(DF))`. That way we have some data to play with right away. — Roman Luštrik, Aug 22 '18 at 05:28
Thanks! I've done that, there are quite a few unrelated variables included there. Let me know if it's okay. — ZoeR, Aug 22 '18 at 05:42
scratch that, I was being lazy, I have just selected the relevant vars :) — ZoeR, Aug 22 '18 at 05:47
Thank you for responding the my hail and improving your question. — Roman Luštrik, Aug 23 '18 at 06:14

score 0 · Answer 1 · answered Aug 22 '18 at 07:05

For doing the interpolation you can do the following:

library(tidyverse)


df <- structure(list(DF.Date = structure(c(1357296180, 1358021760, 
                                           1357991280, 1357629540, 1357287540, 1358158260), class = c("POSIXct", 
                                                                                                      "POSIXt"), tzone = "UTC"), DF.Irradiance = c(1055.64, 0, 975.22, 
                                                                                                                                                   0.25, 953.33, 1017.62), DF.AirMass0 = c(0.0279842300461168, 0.0195227712650674, 
                                                                                                                                                                                           0.0194685363977006, 0.0384509898910604, -0.142855127416022, 0.0322758934562728
                                                                                                                                                   )), row.names = c(NA, 6L), class = "data.frame")
df_new <- approx(x = df$DF.AirMass0, y = df$DF.Irradiance,
                 xout = seq(min(df$DF.AirMass0), max(df$DF.AirMass0), 0.0001)) %>% 
  bind_cols()

This will interpolate your x by 0.0001. After doing this, retrieving the value you are looking for should be straightforward.

Dan · Accepted Answer · 2018-08-22T11:59:32.037

Not the most elegant solution, but this uses dplyr and magrittr. First, I define your data frame.

# Data frame called df
#                  Date Irradiance    AirMass0
# 1 2013-01-04 10:43:00    1055.64  0.02798423
# 2 2013-01-12 20:16:00       0.00  0.01952277
# 3 2013-01-12 11:48:00     975.22  0.01946854
# 4 2013-01-08 07:19:00       0.25  0.03845099
# 5 2013-01-04 08:19:00     953.33 -0.14285513
# 6 2013-01-14 10:11:00    1017.62  0.03227589

Next, I load the relevant libraries.

# Load libraries 
library(dplyr)
library(magrittr)

Here, I create a function that takes an air mass and given x value (i.e., 0.9804), creates a reference data frame (i.e., x equal to 0.5, 0.75, and 1 and corresponding y values), then creates a function that will estimate y based on x through linear interpolation.

# Calculate pressure correction based on air mass and x
pres_cor <- function(m, x){
  # Create reference data frame
  ref_df <- data.frame(x_ref = c(0.5 , 0.75 , 1),
                       y_ref = c(1.68219 - 0.03059 * m + 0.000890 * m^2, 1.248274 - 0.011997 * m + 0.000370 * m^2, 1))

  # Create function for interpolation
  int_fun <- with(ref_df, approxfun(x_ref, y_ref))

  # Return value at given x value
  int_fun(x)
}

Finally, I apply this to each row of your data frame using the pipe operator (%>%), rowwise, and mutate from dplyr, and the compound assignment pipe (%<>%) from magrittr.

# Use function for each row
df %<>% 
  rowwise %>% 
  mutate(y = pres_cor(AirMass0, 0.9804))

This gives the following:

# # A tibble: 6 x 4
#   Date                Irradiance AirMass0     y
#    <dttm>                  <dbl>    <dbl> <dbl>
# 1 2013-01-04 10:43:00    1056.     0.0280  1.02
# 2 2013-01-12 20:16:00       0      0.0195  1.02
# 3 2013-01-12 11:48:00     975.     0.0195  1.02
# 4 2013-01-08 07:19:00       0.25   0.0385  1.02
# 5 2013-01-04 08:19:00     953.    -0.143   1.02
# 6 2013-01-14 10:11:00    1018.     0.0323  1.02

Note that y values look the same due to rounding, but are not upon closer inspection.

# df$y
# [1] 1.019438 1.019446 1.019446 1.019429 1.019600 1.019434

Perfect. Ah, I should have been using approxfun. Thank you very much for your time and help! :) — ZoeR, Aug 22 '18 at 23:54

Having trouble conceptualising how to implement a linear interpolation in R

2 Answers2