We consider a strongly NP-hard Euclidean problem of nding a subsequence in a nite sequence under the criterion of the minimum sum of squared distances from the elements of sought subsequence to its geometric center (centroid). It is assumed that the sought subsequence contains a given number of elements. In addition, sought subsequence has to satisfy the following condition: the difference between the indexes of each previous and next points is bounded with given lower and upper constants. We present an approximation algorithm for the problem and prove that it is a fully polynomial-time approximation scheme when the space dimension is bounded by a constant.
In this paper we study one strongly NP-hard problem of searching a subsequence in a nite sequence of points from Euclidean space. Our aim is to provide an approximation algorithm for this problem.
This work was motivated by poor study of the problem. The problem is also interesting because of its importance for applications, in particular, for mathematical problems of time series analysis, approximation problems and also for applications dealing with data mining problems (see, for example, [1{4] and references therein).
The paper is organized as follows. In the next section the formal de nition of the problem under study is given; an example of application (origin) of the problem is also presented. In Section 3, we provide a review of the previous results and announce the obtained algorithmic result. Basic de nitions and statements that provide necessary elements to prove the properties of the proposed algorithm are presented in Section 4. Finally, in Section 5 we construct an approximation algorithm for solving the problem under study and prove that our algorithm is a fully polynomial-time approximation scheme (FPTAS) when the space dimension is xed.
Copyright ⃝c by the paper's authors. Copying permitted for private and academic purposes. In: A. Kononov et al. (eds.): DOOR 2016, Vladivostok, Russia, published at http://ceur-ws.org
Everywhere below R denotes the set of real numbers, ∥ ∥ denotes the Euclidean norm in Rq.
Problem under consideration has the following formulation (see [5], [6], [7]). Problem 1 (Finding a subsequence in a sequence). Given a sequence Y = (y1; : : : ; yN ) of points from Rq and positive integer numbers Tmin, Tmax and M > 1. Find a subset M = fn1; : : : ; nM g N = f1; : : : ; N g of indexes of the sequence Y elements such that
F (M) = ∑ ∥yj
y(M)∥2 ! min;
j2M
(
Tmin nm nm 1
N; m = 2; : : : ; M; on the elements of the tuple (n1; : : : ; nM ).
Problem 1 has the following interpretation (see [5], [6]). There is a time series containing N measurements y1; : : : ; nN of q numerical characteristics of some objects. Each measurement result in the time series has an error, and no correspondence is known between the elements of the time series and the objects. Some of this objects have identical characteristics (or one can say that in the time series there are several measurements of one signi cant object). Other objects are distinguished and have different characteristics (or one can say that in the time series there are some measurements that are treated as "trash" which could be presented in this time series). The number of measurements for identical objects is known. In addition, it is known that the time interval between every two consequent results of measuring characteristics of the identical objects is bounded from above and below with some constants Tmax and Tmin. The characteristics of identical objects in contrast to the characteristics of other objects have basic information value. It is required to nd the subsequence of measures which corresponds to the identical objects using the criterion of minimum sum of squared distances and to estimate the characteristics of these objects (taking into account the measuring errors in the data). 2
Problem 1 is among poorly studied discrete optimization problems. A special case of
this problem when Tmin = 1 and Tmax = N is equivalent [5] to the strongly NP-hard
problem of searching points subset. In this case at the input there is a set of points
instead of a sequence and time-related constraints (
First, let us recall the results obtained for the searching subset problem because it is a simple particular case of Problem 1.
In general, the searching subset problem is strongly NP-hard [8]. But in the case of xed dimension q of the space this problem could be solved [9] in O(N q+1) time.
Moreover, by this time for the searching subset problem the following algorithms have been presented. In [10] a 2-approximation polynomial algorithm is proposed, it's time complexity is O(qN 2). A polynomial-time approximation scheme (PTAS) is substantiated in [11]. This scheme allows to solve the searching subset problem for arbitrary relative error " in O(qN 2="+1(9=")3=") time. For the case of xed dimension q of the space and integer-valued points coordinates an exact pseudo-polynomial algorithm is constructed [12]. The running time of this algorithm is O(N (M D)q), where D is maximal absolute input point coordinate value. It is established [13], that the searching subset problem has no FPTAS, unless P=NP. In cited paper, an FPTAS was proposed for the case of xed space dimension. This scheme allows to solve the problem with arbitrary relative error " in O(N 2(M=")q) time.
By this time for the considered Problem 1 the following results were obtained. First, we should note that as far as the Problem 1 is generalization of strongly NP-hard subset problem, there are neither exact polynomial-time, nor pseudo-polynomial-time algorithms or FPTAS schemes, unless P=NP.
The case when Tmin and Tmax are parameters of Problem 1 is analyzed in [5]. In this work the authors showed that this problem is strongly NP-hard for any Tmin < Tmax. In the trivial case when Tmin = Tmax this problem can be solved through a polynomial time.
In [6] a 2-approximation polynomial-time algorithm is proposed; the running time of the algorithm is O(N 2(M N + q)). In the case of Problem 1 with integer components of the sequence elements and xed dimension q of the space in [7] an exact pseudopolynomial-time algorithm is substantiated. This algorithm nds an optimal solution of Problem 1 in O(N 3(M D)q) time.
Among the highest interest is the question of approximability of Problem 1. In particular, the question of FPTAS construction for the special case of Problem 1 (subclass of problem) has a great importance because in general case such a scheme does not exist. In the present work such scheme is presented for the case of xed space dimension.
The main result of this work is an approximation algorithm which allows to nd a (1 + "){approximate solution for arbitrary relative error " > 0 in O(N 2(M (Tmax Tmin + 1) + q)(√ 2"q + 1)q) time. If the dimension q of the space is xed then the time complexity of our algorithm is equal to O(M N 3(1=")q=2), and it implements an FPTAS. 3
In order to substantiate our algorithm we'll need a few basic assertions and auxiliary problem with exact polynomial-time algorithm for this problem.
Lemma 1. For an arbitrary point x 2 Rq and a nite set Z Rq it is true that ∑ ∥z z2Z x∥2 = ∑ ∥z z2Z z∥2 + jZj ∥x z∥2; where z is a centroid of set Z. Lemma 2. In assumptions of Lemma 1, if some point u 2 Rq is closest (by distance) to the centroid z of set Z among all points from this set, then ∑ ∥z
Both lemmas are well-known. Their proofs are presented in many publications (for example, in [12], [13]).
Lemma 3. Let
S(M; x) = ∑ ∥yn x∥2; x 2 Rq; M
N ;
(
Sx(M ) = ∑
∥yi
i2M
Finally, combining (
This assertion could be easily veri ed by differentiation of S over x and also follows from Lemma 1.
In addition to Lemma 3 for an arbitrary xed point x 2 Rq the restriction of the function S(M; x) to M we denote as Sx(M) and argument of its minimum we denote x as M .
Lemma 4. Let M be the optimal solution of Problem 1, and y(M ) be the centroid
of the set fyiji 2 M g. Then for any point x 2 Rq the following inequality holds
Proof. Let y(Mx) = jM1xj ∑i2Mx yi be a centroid of fyi j i 2 Mxg. Then from de
nitions (
F (Mx) = ∥yi y(Mx)∥2 ∥yi x∥2 = Sx(Mx): In addition, from the de nition of the set Mx we have Further, Lemma 1 applied up to the point x and the set fyi j i 2 M g implies Sx(M ) = ∥yi x∥2 = ∥yi y(M )∥2 + jM j ∥x y(M )∥2 + jM j ∥x
Lemma 4 shows that quality of feasible solution Mx obtained by some point x 2 Rq could be estimated via distance from this point to (unknown) optimal centroid y(M ). The closer considered point x up to the optimal centroid, the less absolute approximation error of an obtained solution.
Lemma 5. Let M
be the optimal solution of Problem 1, and let be the point of the optimal set fyiji 2 M g closest to its centroid, then t = arg
min
y2fyi j i2M g
∥y
for any i 2 M . Summing up both sides of this inequality for all i 2 M
we obtain
∥t
∥yi
Lemma 6. Assume that conditions of Lemma 5 are held. Then for Mx to be a (1+"){
approximate solution of the Problem 1 for xed " > 0, it is enough to take such point
x that satis es the inequality
∥x
t
F (M )
Proof. Let y(Mt) = 1 ∑i2Mt yi be a centroid of the set fyi j i 2 Mtg. From the
fact that Mt = arg mijnMMtj St(M) and de nitions (
⊔⊓
(
St(Mt)
St(M ): Further, since the set fyi j i 2 M g and the point t satisfy the conditions of Lemma 2, the following estimate holds ∑ i2M ∥yi
t∥2 and, therefore,
St(M ) = ∑ ∥yi
t∥2 i2M 2 ∑
i2M 2 ∑ i2M ∥yi
y(M )∥2
∥yi
y(M )∥2 = 2F (M ):
Combining (
Sx(M) = ∑ ∥yi
x∥2 ! min;
i2M
while elements of the tuple (n1; : : : ; nM ) satisfy the constraints (
A dynamic programming scheme is presented in the next lemma and its corollary. This scheme allows to nd the optimal solution Mx of Problem 2. The presented scheme is based on results from [6], [14] and given here for completeness. Lemma 7. For any positive integer M > 1, such that (M 1)Tmin N 1, and for arbitrary point x 2 Rq the optimum Smxin = minM Sx(M) of Problem 2 could be found as
Smin = min SMx (n);
x
n2!M
(
⊔⊓ (16) where the values of the functions SMx (n), n 2 !M , are calculated using the following recurrence formulas
Smx(n) = { ∥yn ∥yn x∥2; x∥2 + m 1(n) = {j j maxf1 + (m Corollary 1. Elements n1x; : : : ; nxM of the optimal tuple Mx could be found by the formulas: nxM = arg min SMx (n);
n2!M x nm 1 = arg
min n2 m(nxm)
Smx(n); m = M; M 1; : : : ; 2: (18) (19) Let us show the algorithm implementing above scheme in a step-by-step description. A l g o r i t h m A1.
Input : sequence Y, point x, numbers Tmin, Tmax and M .
Step 1. Compute the values ∥yn x∥2 for each n 2 N .
Step 2. Using formulas (17), calculate the values Smx(n) for each n 2 !m while m = 1; : : : ; M .
Step 3. Find the minimum Smxin of the objective function Sx using (16) and the optimal tuple Mx = (n1x; : : : ; nxM ) by formulas (18), (19).
Output : tuple Mx = (n1x; : : : ; nxM ).
In [6], it was proved that the algorithm A1 nds an optimal solution of Problem 2 in O(N (M (Tmax Tmin + 1) + q)) time. In this expression, the value Tmax Tmin + 1 is at most N . Therefore, the running time of the algorithm is estimated as O(N (M N + q)). 4
Proposed approximation algorithm for Problem 1 can be schematically described as follows. For all points from the input sequence we de ne a special box (cube with center in this point) such that at least one of these boxes contains unknown centroid of the optimal subsequence. By the given value of the relative error we construct multidimensional grid with uniform step by all coordinates. For all nodes of the grid we solve the auxiliary problem using dynamic programming scheme, and obtain a feasible solution of Problem 1. Finally, among all obtained feasible solutions we select one with minimal value of the objective function of Problem 1.
For an arbitrary point z 2 Rq and positive numbers h and H we de ne a set of points
D(z; h; H) = fd j d = z + h(j1; : : : ; jq); ji 2 Z; jh jij
H; i = 1; : : : ; qg: Note that points from this set are placed in the nodes of the uniform multidimensional rectangular grid with size 2H and step h between nodes. The center of this grid is at the point z.
For the number of nodes of this grid the following bound holds jD(z; h; H)j (2 Moreover, Lemma 6 de nes the condition on the value of the grid spacing, wherein there is a node closed enough to the centroid of the optimal subsequence (in sense of guaranteed error "). Therefore, de ne the grid spacing as follows:
Theorem 1. For an arbitrary " > 0 Algorithm A of Problem 1 in O(N 2(M (Tmax
Tmin + 1) + q)(√ 2"q + 1)q) time.
nds a (1 + "){approximate solution Proof. Let t = arg miny2fyi j i2M g ∥y y(M )∥ be the point from the set fyi j i 2 M g, that is closest to the centroid of this set. If for this point on step 3 we have equality F (Mt) = 0, then the set Mt is an optimal solution of Problem 1, because the value of the objective function F is always nonnegative.
Now we consider the second case, when F (Mt) > 0. (20) (21)
In accordance with Lemma 5, there is an inequality (
Take the point d = arg
min d2D(t;h;H) ∥ d from the grid, that is closest to the optimal centroid. Since distance from y(M ) to the closest node d from the grid D(t; h; H) does not exceed hpq 2 , we obtain the estimate ∥d Thus, the point d satis es the conditions of Lemma 6. Therefore, the set Md is a (1 + "){approximate solution of Problem 1.
Note, that either on step 3 optimal solution will be found or the point d and the set Md will be considered during operation of the algorithm in step 5. Therefore, at the end the algorithm is guaranteed to obtain at least (1 + "){approximate solution.
Let us now estimate the time complexity of the algorithm. On step 1 to solve
auxiliary problem O(N (M (Tmax Tmin + 1) + q)) operations is required. Step 2 requires
O(qN ) operations, and step 3 could be performed in O(
To generate the grid D(y; h; H) on step 4 it needs to perform O(q jD(y; h; H)j) operations. The cost of constructing the sets Md on step 5 and computing the values F (Md) is equal to (as on step 1) O(N (M (Tmax Tmin + 1) + q)).
As a result, for each of N points from the sequence Y performing steps 1{5 required O(N (M (Tmax Tmin + 1) + q) jD(y; h; H)j) operations. Finally, on step 6 to choose the best solution it needs to perform O(∑y2Y jD(y; h; H)j) operations.
It remains to note that the size of the grid D(y; h; H) is bounded by jD(y; h; H)j (2 1) + q)(√ 2"q + 1)q): Therefore, the total time complexity of the entire algorithm is O(N 2(M (Tmax
⊔⊓
Let us show that if the dimension q of the space is xed, then presented algorithm implements an FPTAS. Indeed, if " 2 (0; 2q], then one of the factors in nal complexity can be estimated as follows (√ 2q " + 1 )q 2q (√ 2q
)q "
q
3q q (
Acknowledgments. This work was supported by the RFBR, projects 15-01-00462, 16-31-00186 and 16-07-00168.