构造字符串的排列,使其跳回最少。

16

寻找一个多项式时间算法或证明以下问题的NP难度:

给定两个字符串 s1=a1, a2,...,aks2=b1,...,bk,其中 s2s1 的随机排列。

现在我们想要用 s2 构建出 s1。构建的过程如下:

s2 中选择一个与 a1 相等的字母并将其删除。

继续进行,选择字母 a2 并将其删除,以此类推,直到 s2 为空。

请注意,同一字母可能会在 s2 中出现多次。设构建的序列为 C = c1, c2,...,ck,使得 C = s1。我们定义 l 为我们需要在 s2 中跳回去选择下一个字母的次数。

例如,如果我们选择了c3c4,但是在原始的s2中,c4的索引小于c3的索引,那么我们将l增加1。
我们的任务是找到最优序列,使得l最小。例如,如果我们有s1=abacs2=acab,输出应该为1,因为我们可以从s2选择第二个"a",然后选择"b",然后跳回并选择第一个"a",最后添加"c"。
我不知道如何在多项式时间内解决这个问题。我想也许有一些方法可以计算出完美匹配并读取最优序列,但我不确定。到目前为止,我只有一个指数级别的算法。
指数级别的算法如下(不确定是否正确,也不知道如何测试它):
def solvenaive(s1, s2, curr_ind):
    if len(s1) == 0:
        return 0
    first_s1 = s1[0]
    
    vorkommen = findOccurrences(s2, first_s1)
    results = []

    for i in vorkommen:
        new_s1 = s1[1:]
        new_s2 = stringPop(s2, i)
        res = solvenaive(new_s1, new_s2, i)
        
        if curr_ind > i:
            results.append(res+1)
        else:
            results.append(res)
    
    return min(results)

1
res 是一个列表,res + 1 不起作用,尽管你的递归算法的一般思路似乎是正确的。至于测试它,是什么阻止你实现 findOcurrences() 和 stringPop() 函数呢? - Swifty
@Swifty,看起来确实有点混乱,但是每个迭代实际上都返回一个数字(min(results)),所以res + 1会起作用。 - codingStarter
1
你的算法已经实现了一半的分支定界。你使用递归正确地对所有可能性进行了枚举。但你从未“剪枝”任何子计算。这正是问题所在。也就是说,你应该估计从一个子调用可以达到的最小“res”,如果那个最小值并不比你已经获得的最佳“res”更小,就立即停止。 - chrslg
顺便提一下,如果您的函数除了最小返回次数之外还返回“路径”,那么可能会更好。 并且,在@ chrslg的评论上进行扩展,也许将 'vorkommen' 分成“ahead_indices”(>= curr_ind)和“backward_indices”(<curr_ind),先处理 ahead_indices,在 res = 0时立即返回0,在最小结果为1的情况下处理完它们后,返回1,然后处理 backward indices,在 res = 0时返回1. 我看不出您还能做什么优化,但也许我漏掉了什么。 - Swifty
从我所了解的情况来看,NP问题必须在多项式时间内可验证。我没有证据,但我的直觉是我们没有独特的方法可以在多项式时间内验证任何一个解决方案,除非对一般的多项式时间解决方案进行正确性证明。这使得这个问题超出了NP范畴,并且根据定义成为NP难题。 - גלעד ברקן
显示剩余2条评论
3个回答

6

为了澄清我在分支限界评论中的观点,这是我的看法。

它显然比你的实现更快(同时保持“天真”的纯Python实现。例如,在我们明确想要“变异”字符串时使用不可变字符串:D)。然而,它仍然显然不是多项式的。只是一个更快的指数。

这在大多数情况下都是分支限界的情况。

import itertools
import random
import numpy as np
import timeit
import matplotlib.pyplot as plt

W=10 # Size of a word

# Generate a random word
def randomWord():
    return ''.join(np.random.choice(list("abcdef"),W))

# And a random shuffling of it
def shuffle(s):
    return ''.join(random.sample(s,len(s)))

# Quick'n'dirty implementation of the function you are using
# It would have been better, btw, if you had provided them, as
# needed to create a minimal reproducible example
def findOccurrences(s, c):
    return [pos for pos, char in enumerate(s) if char == c]

def stringPop(s, i):
    return s[:i]+s[i+1:]

# Your function
def solvenaive(s1, s2, curr_ind):
    if len(s1) == 0:
        return 0
    first_s1 = s1[0]
    
    vorkommen = findOccurrences(s2, first_s1)
    results = []

    for i in vorkommen:
        new_s1 = s1[1:]
        new_s2 = stringPop(s2, i)
        res = solvenaive(new_s1, new_s2, i)
        
        if curr_ind > i:
            results.append(res+1)
        else:
            results.append(res)
    
    return min(results)

# Mine. Almost the same.
# But I add 2 parameters
# bestScore which is what we have seen best
# and sofar, a partial score 
def solvecut(s1, s2, curr_ind, bestScore, sofar=0):
    if len(s1)==0: # If we reach a "leaf" of recursion, then
                   # return the score. (sofar is a real score in that case, not a partial)
        return sofar
    # Otherwise, we will recurse... unless it is a lost cause
    if sofar>=bestScore: # No need to try better, we are already over our best
        return 1e9

    first_s1 = s1[0]
    
    vorkommen = findOccurrences(s2, first_s1)


    for i in vorkommen:
        new_s1 = s1[1:] # I let that one here in the loop because I don't want to optimize anything other that my "cuts"
                        # so that improvement are clearly imputable to those cuts
                        # but, obviously that could have been done outside the loop
        new_s2 = stringPop(s2, i)
        if curr_ind>i:
            res=solvecut(new_s1, new_s2, i, bestScore, sofar+1)
        else:
            res=solvecut(new_s1, new_s2, i, bestScore, sofar)

        if res<bestScore:
            bestScore=res

    return bestScore # Sometimes we'll return that just because we did nothing best, but it doesn't matter


# Test if result are the same on some examples

for i in range(20):
    w=randomWord()
    s=shuffle(w)
    sc1=solvenaive(w,s,0)
    sc2=solvecut(w, s, 0, 1e9)
    print(w, s, sc1, sc2)

# Timeit
# Note that timing depends a lot of the word (letter repetition, etc.). So it is important to do them with the same words
# Especially since this is not very fast, so we need to timit with number=1
W=17
w1=randomWord()
s1=shuffle(w1)
print(timeit.timeit(lambda: solvenaive(w1, s1, 0), number=1))
print(timeit.timeit(lambda: solvecut(w1, s1, 0, 1e9), number=1))

结果

cbaaecffae aaebcffcae 2 2
eedeafdffb ffdbeefaed 3 3
aedefdceba adeefcabde 4 4
aaafccaafb aafacaafcb 2 2
eaffdfafef efdffefaaf 2 2
dedbfffbce bcdffbedfe 3 3
fedfcbaeed bedcfaefde 4 4
ffeebfcfab feaebcffbf 3 3
fcdaddbfbf ffddacbfbd 3 3
deffcdeaea eeafdadfce 4 4
cafeeebcdb beaedfbecc 4 4
beefdaeabd edaefbbdea 3 3
ddccbdbdae cdbdbcdead 3 3
bcbababdef decafbbbba 4 4
dfaeceefea edfcefaaee 2 2
abcdcbbdbf fbbadccbdb 4 4
acbedbaefc cbeacaebdf 3 3
aaffacfbde cfaaeafbfd 3 3
edceddebdc dedcecbded 2 2
ecafabdfea fafaeaedbc 3 3
4.060046362923458
0.05187216284684837

因此,我们可以合理地相信,我的代码始终返回与您的相同的结果。

但是时间几乎快了100倍。

然而(我已经尝试过许多W值,不幸的是,尽管有随机性,它仍然是指数级的):

1 cut=6.040092557668686e-06 
2 cut=8.61193984746933e-06 
3 cut=1.5990110114216805e-05 
4 cut=1.5911879017949104e-05 
5 cut=2.897903323173523e-05 
6 cut=3.5561854019761086e-05 
7 cut=5.8827921748161316e-05 
8 cut=0.00018196692690253258 
9 cut=0.0001191610936075449 
10 cut=0.0002572301309555769 
11 cut=0.0008078841492533684 
12 cut=0.0009584939107298851 
13 cut=0.008062224835157394 
14 cut=0.004587493138387799 
15 cut=0.08773919194936752 
16 cut=0.02209200686775148 
17 cut=0.045466484036296606 
18 cut=0.11209587100893259 
19 cut=0.40787983988411725 
20 cut=2.2789966259151697 
21 cut=1.8927371250465512 
22 cut=1.097903796005994 

不太顺畅...但我们可以清楚地看到,时间日志是“带有很多噪音的线性”,而不是值。因此,它仍然是指数的。

还要注意,通过一些思考,我们应该能够估计出比目前更好的可达分数的最小界限(使用目前的意思是我们仍然希望有一个解决方案,而从未回头。例如,如果curr_ind>0,则我们至少可以通过说如果curr_ind>0,则sofar+1必须小于bestScore,否则就没有尝试的用处)

编辑:我忍不住试了一下。 所以这里有一些轻微的优化。

def lessNaive(s1, s2, curr_ind, bestScore, sofar):
    if len(s1)==0:
        return sofar
    if sofar>=bestScore or (curr_ind>0 and sofar+1>=bestScore): # No need to try better, we are already over our best
        return 1e9

    first_s1 = s1[0]
    new_s1 = s1[1:]

    for i in range(len(s2)):
        if s2[i]!=first_s1: continue
        new_s2 = stringPop(s2, i)
        if curr_ind>i:
            res=lessNaive(new_s1, new_s2, i, bestScore, sofar+1)
        else:
            res=lessNaive(new_s1, new_s2, i, bestScore, sofar)

        if res<bestScore:
            bestScore=res

    return bestScore # Sometimes we'll return that just because we did nothing best, but it doesn't matter

这里没有太多的优化

  1. 我将new_s1的计算放在循环外面
  2. 我将s2位置的迭代与搜索这些位置结合起来
  3. 正如我所说,如果curr_ind>0,那么我们知道我们不能达到小于sofar+1。这可能是最重要的优化,即使只有这个+1,因为它削减了整个递归分支,这些分支即将等于最佳得分。

这次我得到的时间

11.236965486081317
0.4795242650434375
0.08659004187211394

[naive、cut、lessNaive的时间]

显然,这是一个更难的情况(因为它与之前的4.06/0.051相同)。但关键是,即使是像这样的微小优化,例如+1的改变,也可以使计时的数量级提高一个数量级!(然而,它仍然是指数的。一个更小的指数,但是指数的。指数的事情在于两个指数之间的比率也是指数的!因此,即使从一个指数到另一个指数,这真的值得。然而,最小的指数当然仍然不是多项式)

计时更新

在与Kelly的讨论后进行了后续测试。使用6个版本,即已经测试了3个版本(naive、cut、lessnaive),并讨论是否使用lru_cache

行动方式是,我随机生成一个W=22的单词和其洗牌。然后使用相同的单词运行所有6个代码(这就是我之前已经做过的事情。由于计时取决于输入的随机性,因此必须使用相同的输入来进行比较)。与以前的计时有所不同(在lru_cache版本之外),我保留了最小/最大/平均计时比率。

以下是在25个单词后的计时(虽然不多,但使用W = 22和我的相对较慢的计算机已经花了2天时间)

i=25  fn=solveNaive  t= 120.7  μ=1254.1 ρ=100.00 [100.00:100.00]
i=25  fn=solveNaiveC t= 120.5  μ=1254.0 ρ= 99.86 [ 98.89:102.20]
i=25  fn=solveCut    t=  12.0  μ=   9.1 ρ=  9.95 [  0.00:  9.95]
i=25  fn=solveCutC   t=  12.1  μ=   9.1 ρ= 10.01 [  0.00: 10.01]
i=25  fn=lessNaive   t=   4.0  μ=   2.1 ρ=  3.29 [  0.00:  3.29]
i=25  fn=lessNaiveC  t=   4.0  μ=   2.1 ρ=  3.30 [  0.00:  3.30]

C后缀表示带有lru_cache修饰的函数。

所以,我必须说,我没有看到Kelly看到的相同结果。 t是每次迭代的时间(这里只有第25个)。μ是平均时间(所有25个示例的总时间除以25)。然后我计算了一个示例的时间与使用“naive”版本的相同示例的时间之间的比率(因此,“naive”行的比率始终为100%)。ρ是这25个示例的这些比率的平均值。在括号中是这25个示例的该比率的最小值和最大值。

结论是

  • lru_cache 平均来说,对于“naive”版本有一定帮助。但是只是非常微小的帮助,并且对于其他版本没有帮助。
  • 平均时间比例约为30倍,而不是100倍。 100倍是幸运的数字。 但是请注意,这些数字也是最坏情况。原因是有些示例运行时间非常长,这些示例的比例最不利。 然而,它仍然是一个30倍的因素。 请注意,最好的情况是如此之高,我可以测量它(因为我选择仅打印2位小数)。 但这意味着最好的情况是lessNaive的计时比'naive'版本快0.005%以上的场景。 这意味着比例可能会高达x20000,有利于lessNaive。 而且,再次强调,永远不会少于x30,有利于lessNaive

因此,时间确实变化很大。但从未到达 lessNaive 更慢的程度。此外,这似乎几乎不可能,因为它是完全相同的代码,只是停止了一些递归。如果根本没有削减任何递归,那么它可能会非常缓慢:我们将支付测试成本,却没有得到任何回报。但是,该测试的成本比代码的其余部分(pop等)要低得多。即使在所有25个示例中,我们从未遇到过 lessNaive 不至少比正常速度快30倍的情况,理论上可能建立一个示例,其中它略微较慢(一个示例,其中没有递归被行 if sofar>=bestScore or ... 终止),但我认为不可能慢50%。


1
当我运行多次时,时间会极度变化。使用W=22,装饰后的solvenaive有时比lessNaive快100倍,有时又慢2倍。因此,我认为显示单次运行时间是毫无意义且具有误导性的。需要更多的统计数据。 - Kelly Bundy
对于lru_cache,我之前并不知道这个函数。但是如果我理解正确的话,当你使用相同的参数调用函数时它会非常有用。因此优化OP版本可能比我的更有效率(尽管它也可以优化我的版本)。 - chrslg
对于时间问题,我没有遇到过。不同的时间,当然有。我自己也提到过。但是我已经在成千上万次实验中运行了代码。我从来没有看到“lessNaive”变慢。这是在您将装饰添加到solveNaive之前还是之后? - chrslg
1
@KellyBundy 是的,但正确测量时间非常棘手:确保在%timeit迭代之间清除缓存。否则,它只会计算一次,然后在每个后续调用中回忆起该结果... - Pierre D
@PierreD 是的,尽管我删除了 for i in range(20): 循环并且 timeits 有 number=1,所以我没有遇到那个问题。我也不会使用更大的 number,而是会使用多个不同的输入。每次都应该从一个空缓存开始,否则先前的使用可能会使它变得更慢。我的做法是将解决方案函数放入外部函数中,这样解决方案函数每次都会得到一个新的缓存。这样的外部函数还可以隐藏额外的参数,理想情况下,函数应该只接受这两个字符串。 - Kelly Bundy
显示剩余8条评论

4

本次探索的主要发现如下:

  1. 难度与字符串长度n无关,而是与在b中产生a的可能路径数有关(对于给定的n,可能的路径数可以从1到很多)。我们提供了一个闭合形式函数来计算路径数numpaths,它是n_c!(即factorial(n_c))的乘积,其中n_c是字符c在任一字符串中出现的次数。
  2. 因此,这样的潜在解决方案可能非常多。
  3. 在这个集合中,也可能有大量等效的解决方案(成本相同)。
  4. 确定某个路径是否最优(除非其成本为0)是困难的,并且需要工作。
  5. 在分支限界方法中,不同下限之间存在权衡,包括它们的功率(边界有多紧)和计算所需的时间。
  6. 没有一种单一的求解器能够支配所有其他求解器:某些优化在某些设置中表现出色,在其他设置中则表现非常差。即使是@KellyBundy的naive_memoized版本,它以节省时间为代价换取无限空间要求,在平均情况下对于m=n//2表现良好,对于m=isqrt(n)则稍逊一筹,并且在病态情况下可能比其他情况慢1000倍。

所有这些都是提示(但没有证明)问题可能是NP难题。


注意

此前版本的帖子错误地声称已经证明了NP难度。我们曾试图将问题表述为一个有向图,其中最小权重的哈密顿路径对应于我们在这里搜索的解。

首先,正如@KellyBundy和@NiceGuySaysHi正确指出的那样,这只能证明我们可以在指数时间内解决该问题。为了证明NP难度,我们需要做相反的事情:找到一个已知的NP难问题,并展示我们可以将该问题的任何实例转化为我们的问题来解决它,同时只增加输入大小的多项式表达式。

此外,正如@גלעדברקן正确指出的那样,图形过程也是错误的:尽管它可以处理简单的例子,但更大的例子会引入不需要的可能解。

因此,我纠正了错误,并完全删除了图形分析。剩下的是其他一些人发现有趣的分析部分。

概述

一些基本观察: - 距离不对称:dist(a, b) != dist(b, a) - 距离与反转字符串的距离相同:dist(a, b) == dist(r(a), r(b)),其中r(s)是s[::-1] - 这立即淘汰了通常的字符串距离(Levenshtein、Jaro、Jaro-Winkler、Hamming等)
建立一些直觉: - 难度与字符串长度无关 - 什么是复杂性的好衡量?路径数量的表达式,numpaths(b) - 2D表示,poswheels编码 - 很难验证路径是否最小 - 很难开发一个好的下限 - 一些解算器(返回成本和找到的路径) - 时间观察。

基本观察

首先是一个用于生成任意大小问题的工具:

def gen(n, m):
    """Generates a random problem with n chars chosen among m"""
    a = ''.join(np.tile(list(ascii_letters[:m]), n // m + 1)[:n])
    a = shuffle(a)
    b = shuffle(a)
    return a, b

# example
>>> gen(6, 3)
('bcbaca', 'baaccb')

然后,我们还采用了@chrslg的优秀lessNaive()函数,并作出一些小的修改:

HIGHCOST = 1_000_000  # unlikely to ever be reached in our problems
# only for path-returning functions:
HIGHSOL = HIGHCOST, tuple()  # the solution we return when we give up
TERMSOL = (0, tuple())  # leaf solution


def lsn2(a, b, best=HIGHCOST, prev=-1):
    """
    reformulation of `lessNaive()`, eliminating `sofar` and only using `best`
    """
    if not a:
        return 0
    if prev > 0 and best <= 1 or best <= 0:
        return HIGHCOST
    c, a = a[0], a[1:]
    for i, bc in enumerate(b):
        if c != bc:
            continue
        if prev > i:
            res = lsn2(a, stringPop(b, i), best - 1, i) + 1
        else:
            res = lsn2(a, stringPop(b, i), best, i)
        if res < best:
            best = res
    return best

现在,我们检查距离是否对称以及其他一些属性:
df = pd.DataFrame([
    gen(n, (n + 1) // 2)
    for n in range(3, 20)
    for _ in range(100 // n)
], columns=list('ab')).drop_duplicates()

df = df.assign(
    a_b=df.apply(lambda r: lsn2(r.a, r.b), axis=1),
    b_a=df.apply(lambda r: lsn2(r.b, r.a), axis=1),
    ra_rb=df.apply(lambda r: lsn2(r.a[::-1], r.b[::-1]), axis=1),
    a_rb=df.apply(lambda r: lsn2(r.a, r.b[::-1]), axis=1),
    ra_b=df.apply(lambda r: lsn2(r.a[::-1], r.b), axis=1),
).set_index(['a', 'b'])

非对称成本的例子:

>>> df[['a_b', 'b_a']].query('a_b != b_a').reset_index()
                      a                    b  a_b  b_a
0                  bbaa                 abab    1    2
1                  baba                 aabb    2    1
2                 acbba                bacab    1    2
..                  ...                  ...  ...  ...
73  cefibfcjhdahbaedgig  cebfjadidebcihhaggf    7    5
74  hebicjhgfcdfedaagbi  edjcdgcbbfhiaghifea    6    7
75  bifghgecifdaachedbj  ifjfcebadchgeahdgib    6    5

然而:

cols = [k for k in df.columns if k not in {'a', 'b'}]
equalities = [
    f'{i} == {j}' for i, j in permutations(cols, 2)
    if (df[i] == df[j]).all() and i < j
]
>>> equalities
['a_b == ra_rb', 'a_rb == ra_b']

dist(a, b) != dist(b, a), 但是dist(a, b) == dist(reverse(a), reverse(b))

这立即淘汰了大多数常见的字符串距离(Levenshtein、Jaro、Jaro-Winkler、Hamming等),因为它们是对称的。

建立一些直觉

首先要观察的是,找到解决方案的难度并不取决于字符串的大小n。当不存在字符重复时,无论n有多大,都只有一条可能的路径。

相反,主要的困难度量是每个字符位置可用的排列数量。

让我们用可能是最低效的求解器来说明:一个遍历所有可能排列的求解器。请注意,这也可以用来检查另一个求解器的正确性(通过检查其返回的解决方案是否是可能的最小路径之一)。我们还引入了一个(显然的)成本度量,仅使用找到的路径作为输入。

但是首先,让我们介绍一种编码问题的方法,使一些分析变得更容易:

def encode(s):
    """encode a string as a dict: {c: i}"""
    d = {}
    for i, c in enumerate(s):
        d[c] = d.get(c, tuple()) + (i, )
    return d

def encode_as_wheels(a, b):
    """
    Encodes a problem as: pos, wheels.
    For each character at position `i` in `a`, that character can
    be found in `b` at indices `wheels[pos[i]]`.

    Each wheel is always ordered.

    Example:
        >>> encode_as_wheels('aba', 'baa')
        ((1, 0, 1), ((0,), (1, 2)))
    """
    wheels = {}
    for i, c in enumerate(b):
        wheels[c] = wheels.get(c, tuple()) + (i, )
    rw = {c: i for i, c in enumerate(wheels)}
    pos = tuple(rw[c] for c in a)
    wheels = tuple(wheels.values())
    return pos, wheels

有了这个,下面是“by_permutation”求解器:

def compute_cost(seq):
    return sum([i > j for i, j in zip(seq, seq[1:])])


def by_permutation(pos, wheels):
    slots = [v for _, v in sorted(encode(pos).items())]
    for val in product(*[permutations(w) for w in wheels]):
        seq = [(i, j) for ii, jj in zip(slots, val) for i, j in zip(ii, jj)]
        seq = [v for _, v in sorted(seq)]
        yield compute_cost(seq), seq

这里是一个简单的例子:

>>> list(by_permutation(*encode_as_wheels('abacb', 'bacab')))
[(2, [1, 0, 3, 2, 4]),
 (1, [3, 0, 1, 2, 4]),
 (3, [1, 4, 3, 2, 0]),
 (2, [3, 4, 1, 2, 0])]

在这种情况下,有一条成本为1的最优路径([3, 0, 1, 2, 4]),两条成本为2的路径和一条成本为3的路径。
让我们比较两个问题:n=15, m=15n=15, m=4
a, b = 'fdombcglihaeknj', 'clbjmohagfdeink'
n, m = len(a), len(set(list(a)))
>>> n, m
(15, 15)

sol = list(by_permutation(*encode_as_wheels(a, b)))
best = min(sol)

>>> best
(8, [9, 10, 5, 4, 2, 0, 8, 1, 12, 6, 7, 11, 14, 13, 3])

>>> len(sol)
1

相比之下:

a, b = 'dabacbbcadccbad', 'baccddcaadbbcab'
n, m = len(a), len(set(list(a)))
>>> n, m
(15, 4)

sol = list(by_permutation(*encode_as_wheels(a, b)))
best = min(sol)

>>> best
(4, [4, 7, 0, 1, 2, 10, 11, 3, 13, 5, 6, 12, 14, 8, 9])

>>> len(sol)
82944

>>> len([s for s in sol if s[0] == best[0]])
600

难度的度量

如果n不能反映问题的难度,那么是什么?从上面的by_permutation函数,我们看到b中路径的数量为:

其中n_c表示在b中出现字符c的次数。

它只与b(或a,因为每个字符串中的n_c相同)有关。因此我们可以以闭合形式计算出来:

from math import factorial
from collections import Counter

def numpaths(b):
    cnt = Counter(b)
    np = 1  # no pun intended
    for k in cnt.values():
        np *= factorial(k)
    return np

以上是示例:

a, b = 'fdombcglihaeknj', 'clbjmohagfdeink'

>>> numpaths(b)
1

a, b = 'dabacbbcadccbad', 'baccddcaadbbcab'

>>> numpaths(b)
82944

>>> numpaths(a)
82944

请注意,当字符在b中连续时,它们不应该以多种方式尝试。一些求解器(例如lsn3)探索了这一点,但并非所有求解器都这样做。因此,我们认为路径数是问题的基本难度的良好估计,当然,启发式算法(简单或更复杂)可以缩小路径数。
因此,这个数字可能是用作时间测量x轴的更好的度量标准。如果时间允许,我们将根据该数字重新审视下面的时间测量部分。
使用poswheels编码的2D表示
为了理解难度,我们构建了plot_wheels()(稍后发布)。其解释是:每个字符都有自己的颜色,在a[x]b[y]中存在字符的位置上有点。一个求解器(默认为xlsn2)找到的最优路径被绘制为黑线:
a, b = gen(6, 3)
>>> a, b
('cabbac', 'accabb')

plot_wheels(*encode_as_wheels(a, b))

让我们看一下之前的情况,我们找到了600个等效解。我们还指定了解决方案求解器为按排列计算,这样我们就可以看到所有等效的解(所有最小成本的600个):

a, b = 'dabacbbcadccbad', 'baccddcaadbbcab'
plot_wheels(*encode_as_wheels(a, b), sol=by_permutation)

我们目前的主要发现如下:

  1. 由于它们的组合性质,可能会有大量的潜在解决方案。
  2. 在这个集合中,也可能存在大量等效的解决方案(成本相同)。
  3. 除非其成本为0,否则决定路径是否最优是困难的并需要工作。

所有这些都是问题“可能”是NP难的提示(但没有证明)。

难以确定路径是否最小

在我写这篇答案的有趣时光中,我尝试了许多优化/启发式方法来进一步消除分支。我发现很难确定路径是否最优。而且很容易开发一个看似很好的“下界”,可以减少很多工作并显著加快搜索速度,只是后来发现(感谢“perfplot”的“equality_check”),它是错误的,并且存在一些边角情况,导致实际路径的成本低于(错误的)下界规定的成本。

也许这种难度(难以确定路径是否最优)是证明NP难度的一条途径?

求解器

我写了几个求解器。有些是基于 pos, wheels 参数工作的,而另一些则直接处理字符串 a, b。有些具有强大的回溯剪枝 lowerbound 估计,而其他一些则追求简单(它们回溯更多,但通常净效果相当,因为它们每个路径需要做的工作要少得多)。

在这里,我将介绍一个非常快速的求解器(在某些情况下比 lessNaivelsn2 快100倍),并且还会返回找到的路径及其成本:

def xlsn3(a, b, best=HIGHCOST, prev=-1):
    """
    When encountering consecutive letters in b, no need to backtrack
    """
    mincost = int(prev > 0)
    if a == b:
        return mincost, tuple(range(len(a)))
    if best <= mincost:
        return HIGHSOL
    c, a = a[0], a[1:]
    bestseq = tuple()
    i0 = max(prev, 0)
    j = -2
    for i, bc in chain(enumerate(b[i0:], i0), enumerate(b[:i0])):
        if c != bc:
            continue
        if i - j != 1:
            deltacost = int(prev > i)
            cost, seq = xlsn3(a, stringPop(b, i), best - deltacost, i)
            cost += deltacost
            if cost < best:
                best = cost
                bestseq = (i,) + seq
                if best <= mincost:
                    return best, adjust_seq(bestseq)
        j = i
    return best, adjust_seq(bestseq)

def adjust_seq(seq):
    if len(seq) < 2:
        return seq
    y = seq[0]
    return (y,) + tuple(e + (e >= y) for e in seq[1:])

这个问题的特殊之处在于一种优化,看起来可以简单地减少许多组合:当两个或更多字母在b中连续时,尝试它们的反向是没有意义的。只有按顺序使用它们,我们仍然保证找到最佳路径,因为先“消耗”高位置再消耗较低位置不会有任何收益。
让我们看一个例子:
b = 'abccbbbbdddea'
a = shuffle(b)
>>> a, b
('baecbddbabcdb', 'abccbbbbdddea')

fig, axes = plt.subplots(ncols=2)
plot_wheels(*encode_as_wheels(a, b), sol=by_permutation, ax=axes[0])
plot_wheels(*encode_as_wheels(a, b), sol=xlsn3(a, b), ax=axes[1])

计时

我们使用优秀的perfplot来绘制一些算法的性能。在本文中,我只介绍了lsn2xlsn3。基本上,lsn2与@chrlsg的lessNaive函数相同(在性能图中几乎无法区分,因此我在此处将其删除)。lsn3xlsn3的仅返回成本的表亲。我们开发了许多其他函数,这里没有展示,具有各种启发式方法以安全地丢弃子分支。其中一些启发式方法通常是各种复杂度的下限,它们使我们放弃不可能产生成本改进的子分支。权衡是截止值的强度(对于给定问题,我们会遭受多少递归)与额外花费在计算界限本身的时间之间。

我们还添加了@KellyBundy的LRU缓存版本OP的“naive”求解器: naive_memoized。这可以通过简单地记住在解决问题的过程中已经探索过的相同子分支来降低时间要求,带来一定的空间(内存)开销。请注意,在naive_memoized内设置新缓存函数的正确实现,以便重复计时不会仅仅记住以前运行的数据。
关于时间的一个重要问题是,对于给定的 n,复杂度可能会有很大的变化(正如上面所述,因为真正的难点来自可能的排列数)。因此,我们使用一些循环来重复测量 perfplot 的运行时间,并绘制这些运行的平均值。

关于naive_memoized的注意事项:通过以时间换取无限制的内存,这个求解器在m=n//2的情况下表现良好,但对于m=isqrt(n)的情况则不太理想(尽管在半对数图中其斜率仍然更好)。

此外,在某些情况下,它的性能可能比lsn31000倍。例如:

a, b = 'bbaababaababaaabbba', 'abbabaaababaabbbaab'
tmem = %timeit -o naive_memoized(a, b)
# 51.7 ms ± 188 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)

tln3 = %timeit -o lsn3(a, b)
# 49.1 µs ± 41.3 ns per loop (mean ± std. dev. of 7 runs, 10,000 loops each)

>>> tmem.average / tln3.average
1052.6

你能否在基准测试中包含 naive_memoized 吗?我认为它会是最快的。 - Kelly Bundy
我对你的陈述感到困惑,即“找到我们的解决方案等同于在图中找到最小权重哈密顿路径”。例如,在这个问题中,给定acbba,abacb,当我们在创建的路径中的第一个b时,图中有边连接到ab,但我们只允许去另一个b。当我们在创建的路径中的第二个b时,我们再次有比我们允许遍历的边更多。 - גלעד ברקן
b1b2 都与 a0a2 相连,但如果我们的路径到目前为止是 a0 -> b1,则不允许使用边 (b1, a2),因为它与原始字符串不符,必须在路径中有第二个 b。我错过了什么? - גלעד ברקן
我从未提到要返回b1。我说即使有可用的边缘,我们也不被允许从b1前往a2,因为结果路径必须有两个连续的“b”。 - גלעד ברקן
@גלעדברקן:你是对的,“collapsed graph”表述是错误的。我已经完全删除了图形分析部分。 - Pierre D
显示剩余9条评论

2

由于多人发现chrslg的非答案很有用,所以我也会在此发布一个。同时也只是一个速度改进(不是要求的“多项式时间或NP难题”),只需使用functools.lru_cache来记忆化原始解决方案。

使用字符串长度为17的多个随机输入进行基准测试:

result original memoized ratio
4 True  1.6349   0.0034   482
5 True  0.5427   0.0045   120
5 True  2.8031   0.0023  1212
5 True  0.1241   0.0015    82
5 True  0.6741   0.0026   255
4 True  0.3481   0.0018   190
5 True  1.3098   0.0025   529
5 True  1.0044   0.0015   688
5 True  0.2801   0.0012   229
5 True  0.6569   0.0026   257
4 True  0.3757   0.0012   326
4 True  0.1771   0.0051    34
3 True  8.1880   0.0034  2439
4 True  1.3146   0.0049   267
5 True  0.3008   0.0022   139
4 True  2.6785   0.0022  1203
3 True  1.6143   0.0097   167
3 True  1.1957   0.0048   247
4 True 10.2205   0.0037  2735

因此,简���地进行记忆化可以使其速度提高数十到数千倍。从多秒到几毫秒。

对于字符串长度为22(仅考虑记忆化版本,因此忽略零和假值):

result original memoized ratio
3 False  0.0000   0.1776     0
6 False  0.0000   0.0537     0
5 False  0.0000   0.1178     0
4 False  0.0000   0.0360     0
6 False  0.0000   0.0637     0
5 False  0.0000   0.1188     0
4 False  0.0000   0.0929     0
4 False  0.0000   0.2369     0
6 False  0.0000   0.0185     0
5 False  0.0000   0.0393     0
4 False  0.0000   0.2396     0
5 False  0.0000   0.0630     0
5 False  0.0000   0.2074     0
6 False  0.0000   0.0834     0
4 False  0.0000   0.2619     0
5 False  0.0000   0.1232     0
4 False  0.0000   0.3286     0
4 False  0.0000   0.1616     0
4 False  0.0000   0.1822     0
4 False  0.0000   0.1093     0
5 False  0.0000   0.0277     0
4 False  0.0000   0.2417     0
5 False  0.0000   0.0614     0
6 False  0.0000   0.0984     0
6 False  0.0000   0.0072     0
5 False  0.0000   0.1572     0
6 False  0.0000   0.0791     0
5 False  0.0000   0.1145     0
4 False  0.0000   0.0819     0
5 False  0.0000   0.0507     0
6 False  0.0000   0.0620     0
4 False  0.0000   0.2658     0
4 False  0.0000   0.0956     0
5 False  0.0000   0.1937     0
3 False  0.0000   0.0777     0
5 False  0.0000   0.2270     0
4 False  0.0000   0.2641     0
4 False  0.0000   0.0861     0
6 False  0.0000   0.0993     0

chrslg 报告原始版本用时 1254 秒,所以记忆化版本的 ~0.13 秒速度大约快了 ~10000 倍(尽管使用的是不同的计算机)。

当字符串长度为 28 时,我们最终看到超过一秒的时间(我估计原始版本在这样的输入情况下需要 数周 的时间):

result original memoized ratio
5 False  0.0000   1.2373     0
5 False  0.0000   0.1386     0
5 False  0.0000   1.2382     0
5 False  0.0000   1.9010     0
5 False  0.0000   1.2455     0
6 False  0.0000   0.5460     0
6 False  0.0000   1.3230     0
5 False  0.0000   1.0667     0
5 False  0.0000   3.5551     0
5 False  0.0000   3.3685     0
8 False  0.0000   0.2763     0
7 False  0.0000   0.5888     0
5 False  0.0000   0.5947     0
6 False  0.0000   1.3308     0
6 False  0.0000   0.8198     0
6 False  0.0000   1.6232     0
6 False  0.0000   0.3384     0
6 False  0.0000   1.0543     0
5 False  0.0000   1.3115     0
7 False  0.0000   1.6151     0
5 False  0.0000   2.2333     0
6 False  0.0000   1.4533     0
5 False  0.0000   0.8959     0
5 False  0.0000   0.7359     0

代码,大多数来自chrslg(在线尝试!):

import itertools
import random
from time import time
from functools import lru_cache

W=17 # Size of a word

# Generate a random word
def randomWord():
    return ''.join(random.choices("abcdef",k=W))

# And a random shuffling of it
def shuffle(s):
    return ''.join(random.sample(s,len(s)))

# Quick'n'dirty implementation of the function you are using
# It would have been better, btw, if you had provided them, as
# needed to create a minimal reproducible example
def findOccurrences(s, c):
    return [pos for pos, char in enumerate(s) if char == c]

def stringPop(s, i):
    return s[:i]+s[i+1:]

# Your function
def solvenaive(s1, s2, curr_ind):
    if len(s1) == 0:
        return 0
    first_s1 = s1[0]
    
    vorkommen = findOccurrences(s2, first_s1)
    results = []

    for i in vorkommen:
        new_s1 = s1[1:]
        new_s2 = stringPop(s2, i)
        res = solvenaive(new_s1, new_s2, i)
        
        if curr_ind > i:
            results.append(res+1)
        else:
            results.append(res)
    
    return min(results)

def solvenaive_memoized(s1, s2):
  @lru_cache(None)
  def solve(s1, s2, curr_ind):
    if len(s1) == 0:
        return 0
    first_s1 = s1[0]
    
    vorkommen = findOccurrences(s2, first_s1)
    results = []

    for i in vorkommen:
        new_s1 = s1[1:]
        new_s2 = stringPop(s2, i)
        res = solve(new_s1, new_s2, i)
        
        if curr_ind > i:
            results.append(res+1)
        else:
            results.append(res)
    
    return min(results)

  return solve(s1, s2, -1)


print('result original memoized ratio')
t0 = time()
while time() - t0 < 30:
  s1 = randomWord()
  s2 = shuffle(s1)
  t = time()
  r1 = solvenaive(s1, s2, -1)
  t1 = time() - t
  t = time()
  r2 = solvenaive_memoized(s1, s2)
  t2 = time() - t
  print(f'{r1} {r1 == r2} {t1:7.4f} {t2:8.4f} {round(t1 / t2):5}')

顺便提一下,我在我的测量中添加了“naive_memoized”,并在文章末尾加入了一些观察结果。 - Pierre D
@PierreD 是的,我看到了,我收到了所有的编辑通知 :-)。 我在想它在n,2时的斜率。也许它在更大的n时会追上? - Kelly Bundy
是的,在'abbabaaababaabbbaab'中,“困难度数字”(可能路径数)为1,316,818,944,000。 对于n=19,它是最大可能的:9! * 10! - Pierre D
@PierreD 一兆条路径52毫秒不错啊 :-) - Kelly Bundy
...但是49.3微秒更好一些;-) - Pierre D

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接