4

I have used smoothing to create two "functions" fd4 and fd6.

fit6 <- smooth.basis(tid6, zbegfor, fdParobj2) 
fd6 <- fit6$fd 

Plot of the curves

I want to measure the L2 distance between them on the interval [0,1], but I haven't been able to find an appropriate way.

||f − g||_2 = sqrt(int(|f(x)-g(x)|^2,0,1))

The best bet has been this one: How to calculate functional L_2 norm using R, but when I use fd6 instead of f <- function(x) x^2, I get the following message:

"Error in fac - fdmat : non-conformable arrays".

I've spent hours trying to find a solution. Please help me!


Now with reproducible code:

library(fda)

# Smoothing of movement pattern without obstacle rescaled to the interval [0,1]
without <- c(22.5050173512478, 22.5038665040295, 22.5171851824298, 22.5368096190746, 
22.5770229184757, 22.6709727229898, 22.8195669635573, 23.0285400460222, 
23.3240853426905, 23.6895323912605, 24.0905709304813, 24.5674870961964, 
25.129085512519, 25.7433521858875, 26.4096817521118, 27.1338935155912, 
27.906416101033, 28.7207273157549, 29.5431756517467, 30.3697951466496, 
31.2214907341765, 32.0625307132683, 32.8786845916855, 33.671550678219, 
34.4449992914392, 35.1852293010227, 35.8866367048324, 36.5650863548079, 
37.1776116180247, 37.7706354957587, 38.3082855431959, 38.8044130844639, 
39.2471137254193, 39.6193031585418, 39.9685683244076, 40.2345560551869, 
40.4394442661545, 40.5712407258558, 40.6905311089523, 40.712419802203, 
40.6704560575084, 40.5583379372846, 40.3965425630546, 40.1443139907057, 
39.8421899334408, 39.4671160834355, 39.018733225651, 38.5381390971577, 
38.035680135599, 37.4625783280288, 36.8649362406917, 36.2320264206665, 
35.5599736527209, 34.8983871226943, 34.2058073957721, 33.4893682831911, 
32.7568501019309, 32.0241649500974, 31.3036406455137, 30.587636320768, 
29.8962657607091, 29.2297665999702, 28.6003939337949, 28.0003531206639, 
27.433551463149, 26.9088532545635, 26.4265682839796, 25.974193299003, 
25.5553146923473, 25.1701249455904, 24.8107813804098, 24.4776168601955, 
24.167582682288, 23.8726502760669, 23.589703789663, 23.3222235336882, 
23.0616248799115, 22.8185342685607, 22.6767541125512, 22.6567795841271, 
22.6488510112824, 22.6436058079441, 22.6391304188382)

timewithout <- (1:length(without))/length(without) # For scaling
splineBasis = create.bspline.basis(c(0,1), nbasis=25, norder=6) # The basis for smoothing
basis = fdPar(fdobj=splineBasis, Lfdobj=2, lambda=0.00001) 
fitwithout <- smooth.basis(timewithout, without, basis) # Smoothing
fdwithout <- fitwithout$fd 

# Same but movement is over an obstacle
with <- c(22.4731637093167, 22.4655561889073, 22.4853719755102, 22.4989400065304, 
22.5495656349031, 22.666945409755, 22.8368941117498, 23.0846080078369, 
23.4160560011242, 23.8285634914224, 24.2923085321078, 24.8297004047422, 
25.4884540279408, 26.2107053559, 27.0614232848574, 27.9078055119721, 
28.8449720096674, 29.8989669834473, 30.996962022701, 32.1343108758062, 
33.3286403418359, 34.6364870430171, 35.9105342483246, 37.1883582665643, 
38.467212668323, 39.7381525466373, 41.0395064969214, 42.3095531191294, 
43.5708069740233, 44.7881178787717, 45.9965529977777, 47.1643807808923, 
48.284786275036, 49.3593991064962, 50.3863035442644, 51.3535489662494, 
52.2739716491521, 53.1338828493223, 53.9521101656512, 54.7037562884229, 
55.3593092084143, 55.9567618011946, 56.4768579145271, 56.9251919073806, 
57.2971965985674, 57.5937987523734, 57.8158626068961, 57.9554856023804, 
58.009777126789, 57.9863251605612, 57.8932199088797, 57.6988126618694, 
57.4350394069443, 57.1112025796509, 56.7580579506751, 56.2680669960935, 
55.6963799946038, 55.0574070566765, 54.3592140352073, 53.6072275005723, 
52.7876353306759, 51.9172334605074, 50.9879178368431, 49.9953932631072, 
48.9460707853802, 47.8511977258834, 46.6827266395278, 45.4635999409637, 
44.2633368255294, 43.0386729762103, 41.7880095105045, 40.4834298069985, 
39.1610223705633, 37.9241872458281, 36.7158342529737, 35.5408830466013, 
34.4070964101159, 33.307156473109, 32.2514661493348, 31.2475129673168, 
30.2990631096187, 29.4096423238141, 28.590173995037, 27.8437368908309, 
27.17493959411, 26.5779670740351, 26.0377946174036, 25.5731202027558, 
25.1761397934058, 24.8319659155494, 24.5479180062239, 24.2940808334792, 
24.09388897537, 23.934861348149, 23.7999923744404, 23.6877461628934, 
23.5982309560843, 23.5207597985246, 23.4354446383638, 23.3604065265148, 
23.2819126915765, 23.1725048152396, 23.0637455648184, 22.9426779696074, 
22.8079176617495, 22.69360227086, 22.6622165457034, 22.6671302753094, 
22.66828206305, 22.6703162730529, 22.6715781657376)
timewith <- (1:length(with))/length(with)
fitwith <- smooth.basis(timewith, with, basis) # Smoothing
fdwith <- fitwith$fd 

# Plots for understanding
plot(fdwith, col=2) # Smoothed curve for movement over obstacle
plot(fdwithout, col=2, add = TRUE) # Same but no obstacle

# I have to find the L2-distance between these curves
Community
  • 1
  • 1
Solveig
  • 43
  • 4
  • Do the answers to [this SO question](http://stackoverflow.com/questions/24742677/how-to-measure-area-between-2-distribution-curves-in-r-ggplot2/24743607#24743607) help? – eipi10 May 13 '16 at 18:29
  • Thank you for the link. The function "density" used in the recommended answer uses a different smoothing-technique. Since one of my tasks is to look at the consequences of changing parameters in the original one, it's not adequate. The measured distance should be similar, though. – Solveig May 13 '16 at 20:46
  • If you provide a [reproducible example](http://stackoverflow.com/a/5965451/496488) (e.g., some sample data, the `smooth.basis` function and anything else needed to reproduce your problem), then you're more likely to get help. – eipi10 May 13 '16 at 22:26

1 Answers1

2

First, one can take advantage of the possibility to perform arithmetic operations with fd objects: fdwith - fdwithout. Second, maybe there is a better way to extract values from fd objects at specific points, but this also works: predict(newdata = 0.5, fdwith - fdwithout). So,

sqrt(integrate(function(x) predict(newdata = x, fdwith-fdwithout)^2, lower = 0, upper = 1)$val)
# [1] 9.592434
Julius Vainora
  • 47,421
  • 9
  • 90
  • 102
  • Thank you so much! You wouldn't happen to know if there's a way to make the area between the smoothed curves grey in the plot? – Solveig May 15 '16 at 11:37
  • @Solveig, `f1 <- function(x) predict(newdata = x, fdwith); f2 <- function(x) predict(newdata = x, fdwithout); x <- seq(0, 1, length = 100); polygon(x = c(x, rev(x)), y = c(f1(x), rev(f2(x))), col = "grey")` – Julius Vainora May 15 '16 at 12:08
  • Thank you once again Julius :) – Solveig May 15 '16 at 13:22
  • @Solveig, no problem, in the nearest future this is likely to be helpful to me too as I am also working with functional data. – Julius Vainora May 15 '16 at 13:33